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INTRODUCTION 


TRENDS:  PC-Software,  Users  Manual  and  Documentation  for  Robust  Graphical 
Time  Series  Analysis  of  Long  Term  Surface  Water  Quality  Records 

TRENDS  are  a  set  of  three  PC  computer  programs  named  TRX,  TRl  and  TR2,  for  performing  graphical 
time  series  analyses  of  temporal  trends  in  long  term  surface  water  quality  data  series.  The  programs  were 
developed  to  handle  the  unique  features  of  surface  water  quality  sampling  records,  specifically  irregularly 
spaced  series  of  virtually  instantaneous  observations,  that  are  not  readily  amenable  to  conventional  statistical 
methods  developed  for  regular  time  series.  The  programs  generate  laser  printer  graphical  images  that  are 
readily  interpretable  by  individuals  with  minimal  statistical  training.  Hard  disk  report  files  are  also  produced 
for  more  conscientious  analysts  seeking  deejjer  understanding  of  the  particular  data  series  and  the  underlying 
dynamical  processes. 

0.0  HARDWARE  REQUIREMENTS: 

The  current  edition  of  the  TRENDS  programs  is  designed  for  implementation  on  a  386  IBM 
PC  or  PC-clone  with  a  80387  math  coprocessor  operating  under  DOS  3.31  or  higher  and 
connected  to  an  HP  LaserJet  Series  II  laser  printer. 

The  laser  printer  should  have  at  least  4  megabytes  of  extended  memory  in  order  to  hold  the 
graphics  images  produced  by  TRENDS  programs. 

1.0  GRAPHICAL  ANALYSIS: 

Although  the  time  series  models  described  subsequently  are  numerically  founded,  the  TRENDS  softw2ire 
applies  the  results  primarily  to  the  production  of  graphics  that  are  intended  to  provide  users  with  a  visual 
synopsis  of  the  historical  surface  water  quality  record  under  consideration  and  for  inclusion  in  written  reports. 
Programs  also  generate  a  report  file  containing  numerical  summaries  to  cissist  interpretation  of  the  graphics. 
This  first  version  of  TRENDS  programs  does  not  support  the  formal  statistical  tests  mentioned  further  on. 
Integration  formal  testing  with  the  numericcil/graphical  procedures  of  the  current  version  is  expected  in 
future  editions  of  the  softwjire. 

Since  viewers  retain  graphical  images  far  longer  than  written  verbal  and  tabular  accounts,  sound  principles 
of  visual  and  graphical  analysis  as  espoused  by  Chambers  et  al.  (1983),  Tufte  (1983)  and  Clevelamd  (1985) 
received  considerable  attention  in  the  development  of  the  TRENDS  graphics.  The  current  edition  produces 
generally  sound  graphics  most  of  the  time;  however,  in  instances,  the  axis  scaling  routines  of  the  PLOT88 
plotting  software  used  to  generate  the  plots,  yield  unsatisfactory  results  that  cannot  be  overridden.  When 
the  scaling  problem  occurs,  data  are  confined  to  a  relatively  small  portion  of  the  available  graphing  space. 
PLOT88  also  suiTers  from  poor  screen  graphics  and  a  limited  selection  of  output  devices  for  generating  hard 
copies.  In  future,  it  is  expected  that  superior  plotting  software  can  be  found  to  replace  PLOT88. 

One  criticism  of  the  software  by  users  of  the  prototype  versions  was  that  it  can  rapidly  generate  prodigious 
quantities  of  paper  output.  However,  the  alternatives  Jire  exceedingly  time  consuming.  An  early  interactive 
version  suffered  from  poor  screen  graphics  and  the  general  difficulty  that  the  decision  to  stop  the  analysis 
at  some  particular  point  really  cannot  be  made  until  the  analysis  has  run  to  completion.  It  is  the  author's 
practice  to  file  the  complete  output  and  also  to  edit  repwrt  files  before  printing.  Graphics  required  for 
repwrts  can  be  selectively  regenerated  as  required.  Computer  modelling  has  always  demanded  the  generation 
of  copious  hard  copies.  The  results  of  trial  runs  can  always  be  recycled. 
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2.0  REPORT  OUTLINE: 

The  report  is  structured  as  follows: 

-  Introduction 

-  PART  I:  review  of  time  series  analysis  and  filtering  concepts  for  analysis  of  surface  water  quality  data 

-  PART  II:  instructions  for  running  TRENDS  programs  including 

-  preparation  of  input  data  files 

-  preparation  of  program  control  (.CNT)  files  including  selection  of  control  parameters 

-  PART  III:  instructions  for  interpreting  report  files  generated  by  TRENDS  programs 

-  Appendix  A:  Figures  (compiled  separately) 

-  Appendix  B:  Tables  (attached  to  main  report) 

-  Appendix  C:  Example  Report  Files  (attached  to  main  report). 

Users  are  advised  to  read  the  report  thoroughly  at  least  once  prior  to  attempting  to  run  TRENDS  programs. 
The  material  in  PARTs  I  i:  II  does  overlap  somewhat  and  will  likely  be  rearranged  in  future.  Because  it  is 
unavoidable  that  an  instructional  guide  to  graphical  software  include  copious  numbers  of  example  graphics 
(81  figures),  the  figures  that  accompany  the  text  have  been  compiled  separately;  so  that,  readers  may  readily 
flip  through  the  figures  while  reading  the  text. 

3.0  INSTALLATION  &c  GENERAL  USAGE: 

Installation  is  accomplished  by  merely  copying  the  program  diskettes  onto  the  particular  hard  disk  drive 
where  the  programs  are  to  reside.  It  is  recormnended  that  the  individual  programs,  i.e.  TRX.EXE,  TRl.EXE, 
TR2.EXE,  be  installed  in  their  own  separate  directories,  e.g. 

-  D:\TRX 

-  D:\TR1 

-  D:\TR2 

Input  data  files  should  likewise  be  placed  in  separate  directories. 

Programs  are  then  activated  by  typing  the  respective  program  name,  e.g.  TRX  and  hitting  the  enter  key. 
The  program  will  respond  with 

TYPE  CONTROL  FILE  NAME: 

The  user  then  responds  with  the  appropriate  control  file  name,  e.g. 
GRN03.CNT 

say  for  a  Grand  River  nitrate  series.  Note  that  if  the  .CNT  file  resides  in  the  same  directory  as  the  program, 
simply  entering  the  file  name  is  sufficient.  Otherwise,  the  full  path  ncmie  is  required,  e.g. 
F:\GRAND\GRN03.CNT 

To  work  quickly  it  is  recommended  placing  .CNT  files  in  the  same  directory  as  the  progreims,  particularly 
since  .CNT  files  are  specific  by  at  least  one  or  two  lines  to  the  program  being  run. 

As  for  choosing  amongst  the  three  programs,  there  are  no  hard  rules.  A  trial  run  with  progrzon  TRX  (say 
with  2  'runs'  and  5  'iterations',  see  PART  II)  should  indicate  if  seasonality  is  strong  enough  to  warrant  the 
seasonal  adjustment  model.  If  seasonality  is  not  strong,  progrcmi  TRl  should  be  run.  When  there  is  strong 
seasonality  present  in  the  series,  users  may  wish  to  run  program  TR2. 

To  assist  users  in  getting  acquainted  with  the  softwEire,  several  of  the  sample  data  sets  used  in  this  report 
have  been  compiled  on  a  diskette  along  with  the  program  specific  .CNT  files  which  are  stored  in  separate 
directories. 
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1.0  MOE  SURFACE  WATER  QUALITY  MONITORING  DATA: 

Within  the  Ministry  of  Environment's  [MOE'S]  Water  Resources  Branch  [WRB]  long  term  surface  water 
quality  records  are  generated: 

1.  at  the  river  and  stream  sampling  sites  of  the  Provincial  Water  Quality  Monitoring  Network  [PWQMN]. 

2.  for  various  inland  lake  sites  administered  by  Lirraiology  Section. 

Responsibility  for  the  central  PWQMN  data  base  currently  rests  with  Watershed  Management  Section 
[WMS]  and  WMS  staff  are  expected  to  be  primary  users  of  the  software.  In  reality,  beyond  some  river 
monitoring  activities  specifically  administered  by  WMS,  the  PWQMN  is  comprised  mainly  of  sites  operated 
by  the  6  MOE  regions.  Hence,  the  programs  are  intended  for  use  by  MOE  Regioncd  Office  staff.  Since 
throughout  southern  Ontario,  many  of  the  river  basin  Conservation  Authorities  [CA's]  participate  ax;tively 
in  the  PWQMN  program,  there  may  be  many  CA  staff  for  whom  these  programs  should  be  of  assistance. 

Because  of  their  high  quality,  a  long  term  lake  data  set  provided  by  Limnology  Section  served  as  one  of  the 
principal  test  data  sets  employed  for  software  development. 

Other  potential  WRB  data  series  to  which  the  programs  could  be  applied  might  include  long  term: 

1.  Great  Lakes  tc  connecting  channel  water  quality 

2.  drinking  water  intaJce  quality 

3.  municipal  and  industrial  effluent  quality 

4.  groundwater  quality. 

Since  the  programs  were  developed  specifically  for  typiccil  river  and  lake  quality  records,  applications  to 
other  types  of  water  quality  data  cited  above  is  not  guaranteed.  Occeisionally  MOE  staff  have  recourse  to 
data  generated  by  other  agencies,  e.g.  water  intake  records  from  the  private  sector,  for  which  the  software 
may  also  be  suitable. 

2.0  TREND  ANALYSIS  OF  SURFACE  WATER  QUALITY  TIME  SERIES 

Surface  water  quality  sample  data  obtained  by  MOE,  typically  present  discrete  time  series  of  irregu- 
larly spaced  quasi-instantaneous  observations  .  In  addition,  water  quality  time  series  are  frequently 
characterized  by: 

•  short  records 

•  high  natural  variability  (low  signal-to-noise  ratio) 

•  non-normal  sample  populations  with  high  skew  and  inhomogeneous  vciriance 

•  numerous  extreme  values  or  outliers 

•  seasonal  variation 

•  covciriate  effects,  e.g.  flow,  temperature,  etc. 

•  changes  in  measurement  methodology  over  time 

•  censored,  i.e.  <  or  >  data,  according  to  laboratory  detection  limits. 

These  features  are  frequently  in  violation  of  the  underlying  assumptions  on  which  conventional  statistical 
methods  rely.  Through  the  1980's  considerable  progress  has  been  achieved  with  the  development  of  statistical 
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methods  suitable  for  difBcult  to  handle  environmental  data  sets;  however,  most  of  these  specialized  techniques 
are  not  yet  readily  available  in  the  common  commercial  statisticcil  packages. 

From  a  water  quality  perspective,  time  trend  analyses  of  historical  water  quality  sampling  records  are  desired 
to  ascertain  if  the  ambient  levels  of  constituents  of  concern  have  increased  or  decreased  either  as  the  result 
of: 

•  emerging  problems 

•  remedial  actions. 

Time  trend  analysis  is  a  particular  aspect  of  the  broader  discipline  of  time  series  analysis  that  is  concerned 
with  changes  over  the  course  of  time  that  might  be  observed  in  a  series  of  observations.  Temporal  trend  is 
a  somewhat  arbitrary  concept  dependent  on: 

•  the  time  horizon  of  interest  to  the  analyst,  e.g.  short,  medium  or  long  term  shifts 

•  the  time  horizon  of  the  available  data  series 

•  the  time  scale  of  mecisurements. 

Two  fundamental  approaches  have  evolved  for  time  trend  analysis  of  difficult  surface  water  quality  time 
series: 

•  formal  statistical  tests  for  trend 

•  graphically  oriented  numerical  procedures  that  permit  the  visual  assessment  of  data  series  behaviour. 
Hybrids  of  the  two  approaches  have  also  been  developed. 

3.0  FORMAL  STATISTICAL  TREND  TESTS: 

The  principal  formal  statistical  tests  for  trend  in  surface  water  quality  time  series  are  nonparcimetric  rank 
order  statistics,  i.e.  based  on  data  ranks,  that  dejjend  to  a  much  lesser  degree  on  the  usual  statistical 
assumptions,  e.g.  underlying  normal  distributions,  than  the  common  <-test.  These  tests  include  severed 
variations  of  the  seasonal  Mann-Kendall  test  (Hirsch  et  o/.,  1982;  Hirsch  and  Slack,  1984;  van  Belle  and 
Hughes,  1984)  which  is  functionally  analogous  to  determining  Kendall's  nonparametric  rank  correlation 
coefficient  r  between  water  quality  concentrations  and  time.  More  recently  McLeod  and  Hipel  (1990)  have 
developed  a  test  based  on  Spearman's  nonparametric  partial  p  statistic  that  performs  at  least  as  well  as  the 
SMK  tests.  While  superior  to  paremietric  tests,  e.g.  <-test,  experience  shows  that  these  tests  are  not  perfect 
and  that  several  alternative  procedures  should  be  applied  jointly. 

It  is  importcint  to  note  that  formal  statistical  tests  are  designed  to  assess  the  net  or  global  monotonic  trend 
over  the  time  horizon  of  the  data  series  under  investigation  (note  that  monotone  trend  is  always  in  the 
same  direction  or  equivalently  there  are  no  reversals  of  direction).  While  these  tests  are  useful  for  assessing 
net  long  term  change  in  a  data  record,  they  will  not  reveal  short  term,  e.g.  year-to-ye£ir  variations,  that  in 
many  situations  may  be  of  equal  or  greater  concern  to  the  analyst.  The  formal  test  statistics  also  are  also 
rather  abstract  and  tend  to  deter  analysts  with  limited  statistical  training.  For  these  reasons,  the  graphically 
oriented  numerical  methods  of  the  TRENDS  software  were  developed. 

4.0  TIME  SERIES  MODELS 

To  develop  numerical  time  series  models,  it  is  necessary  to  formulate  a  conceptual  model  of  the  essential 
features  of  the  data.  Let  Z,  represent  be  a  series  of  water  quality  concentrations  obtained  at  chronologic  times 
ti.  Consider  A',  to  be  a  transformed  series,  e.g.  Zi  =  logjo  A',  where  the  family  of  permissible  transforms 
includes  the  null  transform,  i.e.  Z,  =  A',.  In  the  discussion  below,  the  various  model  components  are 
developed  in  an  algebraic  sense.  Details  of  fitting  procedures  are  outlined  in  the  following  section.  For  the 
present,  visualize  both  the  trend  and  seasonal  components  described  below  as  smooth  lines  drawn  through 
noisy  data  clouds  as  in  Figures  la  and  lb.  Fitting  procedures  and  graphical  features  are  explained  in  the 
later  sections. 

4.1  Simple  Model: 

Progrcim  TRl  performs  simple  trend  analysis. 
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For  a  series  of  noisy  surface  water  quality  measurements  A',  obtained  at  times  t,,  the  simplest  conceptual 
model  for  the  purpose  of  trend  analysis  is: 

A-.  =  T.  +  f.  (1) 

where  T,  is  the  trend  component  and  £,  is  an  irregular  residual  noise  term.  Usually  in  time  series  analysis, 
'trend'  T,  is  considered  to  be  some  measure  of  central  tendency  and  it  is  the  changes  in  T,  over  time  that 
concern  us. 

4.2  Seasonal  Adjustment  Model: 

Programs  TRX  and  TR2  perform  seasonally  adjusted  trend  analysis. 

Many  water  quality  variables  are  known  to  exhibit  regular  annual  periodic  variations  or  seasonal  cycles  in 
which  case  it  is  necessary  to  consider  the  more  complex  seasonal  adjustment  model: 

A'.  =  T.  +  5.  +  £i  (2) 

which  includes  the  seasonal  term  5,.  Models  of  this  general  conceptual  form  are  used  to  produce  the  myricid 
of  'seasonally  adjusted'  economics,  unemployment,  housing  starts  and  other  statistics  regularly  reported  in 
the  news  media. 

From  the  perspective  of  time  trend  analysis,  the  seasonal  component  is  a  nuisance  factor  that  should  be 
eliminated  in  order  to  obtain  reliable  results.  From  a  broader  scientific  perspective,  information  about 
seasonal  variation  can  provide  considerable  insight  into  the  underlying  dynemiic  processes  governing  the 
behaviour  of  water  quality  variables,  and  thus  lead  to  better  conceptual  models  of  water  quality  processes 
and  also  better  monitoring  and  survey  design. 

Suppose  for  a  moment  that  the  trend  term  71  was  constant  over  the  horizon  of  the  data  series,  say  at  the 
mean  level  of  the  series,  i.e.  all  T,  =  T  -  ^"  X,/n.  If  strong  seasonality,  say  a  sharp  rise  over  the  summer 
months,  were  present  in  the  data  series,  we  should  be  able  to  observe  its  presence  by  reordering  the  residucds 
Xi  —  T,  or,  since  Ti  are  constant,  observations  A',  into  the  sequence  A'j  £iccording  to  the  'seasonal  times'  r,- 
at  which  the  observations  occur  within  the  year.  As  before  we  would  consider  defining  the  seasonal  Sj  as 
some  measure  of  the  central  tendency  of  observations  Xj  occurring  in  the  vicinity  of  seasonal  time  r,.  In 
this  respect  the  definition  of  Sj  amounts  to  determining  seasonal  trend  within  the  seasonally  reordered  data 
according  to  the  model 

Xj  =  Sj  +  f,  (3) 

in  direct  analogy  to  model  (1).  Finally  we  reorder  S,  to  chronologic  times  t,  to  obtain  the  seeisonal  comp)onent 
5,  for  use  with  model  (2). 

Note  that  the  seasonal  S,  as  defined  above  is  essentially  a  statistical  measure  of  net  systematic  annually 
recurrent  effects  embedded  in  the  data  series.  In  individual  years  the  seasonal  effect,  e.g.  an  annual  peak, 
may  arrive  earlier  or  later  than  the  net  peak  observed  in  a  record.  Likewise,  in  any  given  year,  the  magnitude 
of  that  peak  may  be  greater  or  lesser  than  the  expected  net  behaviour.  Also  without  further  knowledge  or 
additional  measurements  we  cannot  immediately  associate  Si  with  specific  physiccd,  chemical,  biological  or 
anthropogenic  phenomena. 

4.2.1  Model  Fitting  by  Iterative  Decomposition: 

If  trend  term  T,  is  not  constant,  the  procedure  for  separating  observations  A',  into  Ti,S,  and  f,  rather 
more  difficult  than  for  the  simple  model  (1).  Typically  an  iterative  decomposition  scheme  is  applied  as 
follows: 

1.  with  simple  model  (1),  develop  a  trial  estimate  of  trend  7^' 

2.  obtain  de-trended  residuals  i/,  =  A',  —  T^' 

3.  reorder  de-trended  residuals  i/,  by  time  within  the  year  to  obtain  the  series  i/j 
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4.  from  Uj  and  conceptual  model  (3),  obtain  trial  seasonal  Sj  noting  that  Sj  define  a  relative  seasonal 
comjwnent  that  is  centred  about  0 

5.  obtain  S/  by  reordering  Sj  to  chronologic  times  t,- 

6.  obtain  de-seasonalized  residuals  r;,  =  X,  —  S} 

7.  obtain  a  new  trend  estimate  T^,  by  applying  model  (1)  to  the  de-seasonalized  residuals  r;, 

8.  repeat  steps  2-7  until  estimates  of  T,  and  S,  stabilize. 

Because  the  initial  guesses  at  the  trend  T,^  and  the  seasonal  S/  may  be  rather  crude,  the  iterative  approach 
is  advised.  The  matter  is  discussed  more  thoroughly  in  the  section  describing  program  control  parameters. 
McLeod  et  al.  applied  an  iterative  decomposition  scheme  for  the  analysis  of  monthly  water  quality  time 
series.  The  iterative  smoothing  procedures  developed  for  the  TRENDS  programs  are  also  a  form  of  back 
fitting  (Hastie  and  Tibshirani,  1990). 

Once  the  components  T^^  and  Sj  have  been  successfully  estimated,  the  same  notions  can  be  employed 
to  recompose  the  series  in  alternative  configurations.  For  example,  one  may  adjust  the  entire  series  to  a 
particular  base  level,  e.g.  T,  and  explore  the  seasonal  probabilities  of  exceeding  certain  water  quality  criteria, 
as  per  Bodo's  (1991)  examination  of  seasonal  atrazine  exceedences  of  Canadian  Water  Quality  Guidelines. 
Such  features  will  likely  become  available  in  future  editions  of  the  software. 

N.B.  When  the  logjo  transform  has  been  invoked,  to  be  interpretable  in  the  natural  domain  of  measure- 
ment, the  seasonal  function  must  be  added  to  some  chosen  trend  level  in  the  transform  domain  before  back 
transformation.  This  feature  is  not  currently  available,  will  likely  appear  in  future  editions. 

4.3  Covariate  &:  Seasonal  Adjustment  Model: 

The  following  models  provide  an  appreciation  of  the  limitations  of  current  TRENDS  software 
and  of  what  future  software  development  is  required. 

Although  models  (1-3)  are  the  basis  of  the  TRENDS  software,  it  is  important  to  consider  this  class  of  models 
because  the  successful  determination  of  trend  T,  against  the  high  noise  background  typical  of  many  water 
quality  data  series  involves  eliminating  determinable  systematic  sources  of  variation  to  the  greatest  extent 
possible.  For  many  water  quality  variables,  variation  is  frequently  induced  by  certain  readily  measurable 
covariables.  For  exEimple,  in  rivers  the  instantaneous  flow  rate  exerts  considerable  influence  over  sediment 
transport  capacity;  hence,  variations  sediment  associated  variables  such  as  phosphorus  are  often  largely 
attributable  to  flow  effects.  Likewise  water  temperature  is  a  primary  determinant  of  dissolved  oxygen  levels. 

The  availabihty  of  a  single  coincidentally  observed  covariate  series,  say  flows  Qi,  leads  to  the  possibility  of 
models  such  as: 

A',  =  T,  +  C,  +  S,  +  £,  (4a) 

Xi  =  Ti  +  fi{C,S)  +  £,  (46) 

where  C,  represents  a  systematic  covariate  effect  in  model  (4a)  and  fi(C,S)  in  model  (4b)  represents  a  joint 
covariate-seasonal  effect  that  is  not  necessarily  linearly  additive  as  in  (4a).  Model  (4a)  is  usually  applied  in 
a  manner  analogous  to  the  seasonal  decomposition  scheme  cited  above,  but  without  iterations,  i.e. 

1.  a  concentration-covariate  relationship  is  fitted 

2.  and  then  used  to  determine  covariate  effect  d  as  the  expected  concentration  for  a  given  flow  Qi 

3.  from  which  covariate-adjusted  residuals  A';  —  C,  are  determined, 

4.  then  any  remaining  seasonal  variation  S,  is  determined  on  residual  series  A'i  —  d 

5.  after  which  formal  trend  testing  is  performed  on  the  final  residuals  A',  -  C,  -  S,. 
Work  with  river  quality  series  (McLeod  and  Hipel,  1990;  Hirsch  et  ai,  1991)  has  shown  that: 

•  using  flow  as  a  covariate  can  significantly  reduce  variation  in  the  concentration  series, 
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•  to  varying  degrees,  the  seasonal  variation  for  certain  cases  is  largely  attributable  to  seasonsil  variation 
in  flow, 

•  in  many  other  cases,  after  removal  of  flow  variation,  significant  seasonal  variation  not  attributable  to 
flow  remains  embedded  in  the  data  series. 

At  this  point,  further  work  is  necessary  to  establish  the  best  approach  to  handling  jointly  covariate  and 
seasonal  effects.  Work  with  PWQMN  river  quality  series  suggests  that  for  certain  water  quality  variables 
concentration-flow  relationships  vary  continually  throughout  the  year,  i.e.  seasonally.  Research  is  needed  to 
ascertain  if  the  simpler  linear  additive  model  (4a)  is  sufficient  or  if  a  more  complex  expression  of  concentration 
in  terms  of  a  joint  covariate-seasonal  function  /,(C,  S)  as  in  model  (4b)  is  appropriate.  Also,  it  would  seem 
that  some  form  of  iterative  decomposition  scheme  would  lead  to  a  better  fit  of  either  form  (4a)  or  (4b). 

In  the  future,  once  available,  this  class  of  models  should  be  considered  as  they  embrace  the  elementary  forms 
(1-3)  and  can  be  reduced  as  data  dictate.  A  major  area  of  application  is  the  study  of  river  quality  series  that 
are  highly  reflective  of  non-point  source  (NPS)  pollution  which  is  driven  by  hydrometeorologic  processes  that 
are  highly  random  in  character.  The  random  or  stochastic  nature  of  hydrologic  processes  induces  substantial 
variation  in  the  measured  water  quality  data  series  which  renders  the  detection  of  water  quadity  trends  due 
to  NPS  remediation  measures  a  major  challenge.  A  recent  MOE  sponsored  study  (McLeod  and  Hipel,  1990) 
has  shown  that  very  subtle  trends  can  indeed  be  detected  by  a  combination  of  formed  statistical  tests  and 
flow-seasonal  adjustment  techniques.  Further  work  is  required  to  develop  these  results  into  a  convenient 
softwcire  package. 


5.0  FILTERING  OPERATIONS: 

If  we  consider  model  (1),  practically  the  process  of  defining  Ti  amounts  to  a  smoothing  or  filtering  operation  to 
sepcirate  observations  AT,  into  trend  estimates  T,  and  noise  terms  £, .  Analogously,  in  the  seasonal  adjustment 
model,  the  terms  T;  and  5,  are  determined  by  filtering  operations  performed  on  residual  series  v,  and  r], 
respectively.  Given  the  difficult  characteristics  of  typical  surface  water  quality  time  series,  the  smoothing 
filters  employed  in  the  TRENDS  software  were  selected  to  be: 

•  suitable  for  irregularly  spciced  observations 

•  resistant  to  high-leverage  data,  i.e.  outliers,  such  that  extreme  behaviour  is  absorbed  into  the  noise 
terms  c, 

•  robust  against  asymmetric  data,  i.e.  non-normcd,  data  distributions 

•  independent  of  any  predetermined  trend  model,  e.g.  linear,  step,  monotonic,  etc. 

•  tunable  to  permit  the  examination  of  short,  medium  or  long  term  trend  according  to  the  degree  of 
smoothing  applied. 

The  fundamental  forms  of  filter  employed  in  TRENDS  programs  are  considered  below. 

5.1  Annual  Medians: 

The  annual  median  filter  is  used  by  the  program  TRX  to  determine  the  trend  component  of  the  seasoned 
adjustment  model. 

Although  not  usually  considered  as  such,  a  crude  but  effective  and  readily  implementable  filtering  operation 
that  is  robust  against  outliers  and  distribution  shape,  is  to  determine  a  constant  robust  central  tendency 
measure  such  as  the  median  of  observations  available  for  each  calendar  year.  This  is  known  as  a  'bin  smoother' 
(Hastie  and  Tibshirani,  1990).  The  result  is  a  discrete,  piecewise  smooth  trend  T,  with  discontinuous  step 
jumps  at  the  beginning  (or  equivalently  the  end)  of  each  year  of  the  data  series.  An  example  is  shown  in 
Figure  2  where  for  plotting  convenience  all  sample  dates  have  been  converted  to  decimal  format,  e.g.  April 
24  1985  =  85.311.  This  'annual  median  filter'  has  the  advantage  of  being  very  quick  to  calculate.  To  its 
detriment,  the  median  is  increasingly  unstable  as  measure  of  central  tendency  for  small  (<  10)  sample  sizes 
with  the  consequence  that  the  trend  indicated  for  poorly  sampled  calendar  years  may  be  somewhat  distorted 
or  exaggerated.  On  Figure  2,  this  is  true  for  1981  where  the  indicated  level,  based  on  onlv  4  observations. 
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ranks  amongst  the  3  highest  years.  Analysis  with  the  quasi-continuous  smoothers  discussed  below  support 
the  impression  that  1981  zinc  levels  may  have  been  higher  that  those  of  the  surrounding  years;  however,  the 
suggested  1981  level  is  somewhat  lower  than  that  of  Figure  2. 

5.2  Robust  Smoothing  Filters: 

Programs  TRl  and  TR2  offer  an  array  of  robust  smoothing  filters  derived  from  a  'running  median'  smoother 
developed  for  regularly  spaced  time  series  by  Tukey  (1977)  and  a  "running  line"  x-y  scatterplot  smoother 
develof)ed  by  Cleveland  (1979).  Both  the  Tukey  and  Cleveland  filters  may  be  regcirded  as  'quasi-continuous' 
as  on  the  graphics  to  be  presented  the  filter  output  appears  generally  as  a  smooth  line  traced  through  a 
noisy  data  cloud.  In  reality,  the  methods  yield  a  smoothed  estimate,  e.g.  trend  T,  for  concentration  A',  at 
time  t,,  and  it  is  the  plotting  procedure  of  connecting  the  points  T,  which  produces  the  visual  effect  of  a 
continuous  trend  line.  The  list  of  filters  currently  available  in  programs  TRl  and  TR2  are  given  in  Table 
1.  The  present  discussion  focuses  on  the  concepts  underlying  the  'primitive'  filters  from  which  compound  or 
hybrid  filters  were  developed. 

5.2.1  Tukey  Filters: 

Tukey  (1977)  developed  a  series  of  robust  filtering  procedures  for  tracing  smooth  lines  through  evenly  spaced 
data  series.  The  filters  are  constructed  from  primitive  operations  including: 

-  banning,  a  simple  weighted  average  of  three,  i.e.  (A',  -I-  2A',+i  -I-  A'i+2)/4 

-  repeated  running  medians  of  three,  i.e.  A",  is  replaced  by  the  median  of  (A',_i,  A',,  A'.+i) 

-  splitting,  a  procedure  for  rounding  flat  peaks  and  valleys  comprising  2  or  more  successive  identical 
values. 

The  Tukey  filters  are  designed  for  evenly  spaced  series.  The  degree  of  smoothing  is  not  especially  powerful 
nor  are  these  filters  particularly  well  suited  for  the  irregularly  spaced  water  quality  data.  Nevertheless,  the 
Tukey  filter  (Filter  8),  specifically  an  3RSSH  filter  preceded  by  a  hann,  i.e.  H3RSSH,  is  useful  for  exploratory 
work  as  in  Figure  3  where  the  recurrent  spikey  behaviour  suggests  there  may  be  a  strong  seasonal  component 
in  the  Arsenic  series.  Filter  8  is  cilso  useful  as  a  quick  filter  for  polishing  the  rough  corners  of  the  output  of 
other  filters. 

Tukey's  notion  of  repetitive  filtering,  i.e.  filtering  the  output  of  a  previous  filtering  operation,  is  an  importaint 
concept  exploited  in  the  development  of  higher  order  filters  available  with  programs  TRl  amd  TR2.  The 
concept  of  moving  medians  was  also  exploited  in  the  development  of  the  seasonal  smoothing  filter  described 
below. 

5.2.2  Cleveland's  LOWESS  Filter: 

Cleveland  (1979)  developed  an  algorithm  known  as  LOWESS  to  perform  robust  locally-weighted  regression 
for  the  purpose  of  exploring  nonlinear  relationships  in  x-y  scatterplots.  LOWESS  employs  the  standcird 
linear  least  square  regression  model  y  =  l3o  +  l3ix  +  e,  locally  on  subsets  of  m  <  n  nearest  neighbours  of 
each  of  the  series'  n  data  points  (A',,t,)  where  t,  is  the  chronologic  date.  Moreover,  the  m  local  data  points 
are  weighted  according  to  their  distaiice  from  (A',, t,).  After  a  locally-weighted  regression  estimate  A,  is 
obtained  for  each  of  the  (A",, t,),i  =  l,n  data  points  in  the  set,  additional  robustness  weights  are  defined 
according  to  the  absolute  magnitude  of  the  locally-  weighted  regression  residuals  e,  =  lA",  —  Ai|.  These 
robustness  weights  are  applied  on  a  following  iteration  to  minimize  the  influence  of  outliers.  Note  that  it 
is  the  estimate  A,  that  is  ultimately  used  as  the  trend  term  T,  at  date  i,.  Standard  LOWESS  involves  2 
robustness  iterations.  Bodo  (1989)  and  others  (McLeod  and  Hipel,  1990;  Hirsch  et  ai,  1991)  have  adapted 
the  LOWESS  algorithm  for  tracing  smooth  trend  lines  through  noisy  time  series  of  water  quality  data. 

As  per  Table  1,  in  the  current  edition  of  TRl  and  TR2,  standard  LOWESS  with  2  robustness  iterations  may 
be  performed  with  Filter  1,  while  Filter  9  offers  elementary  LOWESS,  i.e.  locally  weighted  regression  without 
robustness  iterations.  Note  that  LOWESS  with  2  robustness  iterations  amounts  to  3  filtering  operations  each 
with  a  locally  weighted  regression  estimate  developed  at  each  point  of  the  n  points  of  the  data  series.  The 
LOWTSS  algorithm  is  a  computer  intensive  procedure. 
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5.2.2.1  LOWESS  Control  Parameters: 

A  critical  problem  in  applying  LOWESS  to  smoothing  irregular  time  series  is  how  to  develop  systematic 
rules  for  selecting  the  number  of  nearest  neighbours  m  in  a  manner  that  associates  the  choice  with  short, 
medium  or  long  term  trend  definition.  The  solution  developed  for  the  TRENDS  programs  TRl  and  TR2  is 
bcised  requires  two  control  parameters: 

1.  a  typical  annual  sampling  rate 

2.  a  'unit  smoothing  interval'  or  USI. 

The  typical  annual  sampling  rate  is  multiplied  by  the  USI  to  obtain  the  number  of  nearest  neighbours  m 
used  in  LOWESS  and  derivative  smoothing  filters. 

5.2.2.1.1  Typical  Annual  Sampling  Rate  : 

Users  do  not  actually  determine  a  typical  annual  sampling  rate,  but  rather  must  choose  from  two  possible 
ways  of  determining  the  typical  annual  sampling  rate: 

a.  median  rate  =  1  /  median{A,} 

b.  mean  rate  =  1  /  mean{A,} 

where  A,  the  sampling  interval  series  defined  as  A,  =  <,+  i  —  t,:i  =  l,n  —  I.  The  median  sampling  interval 
is  generally  lower  than  the  mean  interval:  hence,  the  median  rate  is  usually  higher  than  the  mean  rate. 
Experience  with  PWQMX  river  quality  series  reveals  that  for  the  vast  majority  of  approximately  monthly 
sampling  sites,  the  median  rate  provides  an  estimate  of  typical  annual  sampling  rate  that  is  unaffected  by  a 
few  long  gaps  in  the  record  and  is  generally  close  to  the  specified  monthly  sampling  objective. 

The  actual  typical  sampling  rate  employed  by  the  programs  TRl  and  TR2  is  determined  by  truncating  the 
estimated  rate  and  adding  1,  e.g.  if  the  calculated  rate  is  10.5  samples,  the  rate  used  by  TRl  and  TR2  is  10 
+  1  =  11,  in  part  for  the  reason  cited  below. 

There  is  cilso  a  technical  reason  for  using  a  slightly  higher  sampling  rate  than  that  which  obtains  most 
years.  Suppose  that  the  median  rate  is  determined  to  be  13  samples/annum  which  is  identically  the  number 
of  nearest  neighbours  m  employed  in  the  LOWESS  algorithm.  Note  that  m  includes  the  point  at  which 
LOWESS  is  being  performed,  i.e.  (X,,t,)  where  ti  is  the  chronological  date.  Suppose  that  (X,,t,)  is  the 
midpoint,  i.e.  the  7th  of  the  nearest  neighbours;  such  that,  there  are  6  points  to  the  left  and  6  points  to  the 
right.  Further  suppose  that  the  data  are  spaced  at  equal  intervals.  Under  this  scenario,  the  two  most  remote 
points  (A',_6,  ti-e)  and  (A'.+e.'i+e)  are  equidistant  from  (A',,  t,)  and  will  be  assigned  local  regression  weights 
of  0,  implying  that  the  local  regression  really  involves  only  the  central  11  nearest  neighbours.  Because  the 
outermost  nearest  neighbours  are  assigned  0  or  very  low  local  regression  weights,  when  performing  LOWESS 
smoothing  with  small  m  as  will  frequently  occur  with  'monthly"  PWQMN  river  records  when  year-to-year 
changes  are  emphasized,  it  is  advantageous  to  have  the  few  additional  newest  neighbours  that  accrue  from 
using  the  median  sampling  rate.  For  example,  with  'monthly'  river  series  the  mean  rate  is  10-11  samples/year 
while  the  medicin  rate  is  12-13  seimples/yecir.  As  a  matter  of  interest,  for  the  Arsenic  series  shown  on  Figure 
3,  the  mean  rate  was  17.2  2ind  the  median  rate  was  19.2.  At  these  levels,  the  effect  in  the  trend  term  estimates 
is  small. 

The  mean  rate  has  been  found  to  be  the  most  appropriate  measure  of  typical  annual  sampling  rate  for  the 
type  of  monitoring  conducted  in  certain  lake  monitoring  programs  for  which  samples  are  obtained  only  from 
May-October  and  a  5  month  record  gap  from  November- April  is  a  regular  feature.  For  these  data  sets,  the 
median  rate  reflects  the  much  shorter  intervals  between  the  May-October  samples  and  may  overestimate  the 
typiccil  cinnual  rate  by  a  factor  or  2-3  or  more.  The  analysis  of  one  lake  record  with  this  type  of  monitoring 
strategy  (see  PART  II:  5.4)  showed  the  mean  rate  was  16.1  samples/year  while  the  median  rate  was  45.6 
samples/year.  On  calendar  year  basis,  over  the  14  years  the  annual  number  of  observations  ranged  from 
6-27  with  a  medicin  15  and  a  mean  15.4. 

The  mean  rate  may  also  be  the  best  choice  for  certain  river  monitoring  records  which  regularly  contain 
winter  record  gaps  due  to  periods  of  ice  cover.  Sampling  rate  information  is  printed  as  a  feature  of  the  input 
routines  of  programs  TRl  and  TR2.  From  a  trial  run,  the  best  choice  of  of  typical  annual  sampling  rate  can 
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be  determined.  If  the  median  rate  exceeds  the  mean  rate  by  a  factor  of  2  or  more  the  mean  rate  is  the  best 
choice.  Also  when  the  ratio  of  median  rate  to  mean  rate  approaches  1.5,  the  mean  rate  is  more  likley  to  be 
the  best  choice.  Diagnostic  programs  are  under  development  to  assist  users  with  these  choices. 

5.2.2.1.2  Unit  Smoothing  Interval: 

The  degree  of  smoothing,  hence  the  focus  on  short,  medium  or  long  term  trend,  is  further  governed  by 
the  choice  of  a  parameter  referred  to  as  the  unit  smoothing  interval  (USI).  To  emphasize  yeaj-to-year 
changes,  the  USI  should  be  set  to  1  year.  For  medium  term  changes,  and  greater  smoothing,  USI  could  be 
set  to  2  or  3  years.  It  is  also  possible  to  set  the  USI  to  fractions,  e.g.  USI  =  .33  year,  which  can  be  useful 
for  exploratory  work  with  program  TRl.  Generally  program  TR2  is  designed  for  examining  trends  with  USI 
taken  as  integral  multiples  of  1  year. 

5.2.3  Modified  LOWESS  (LOWESP)  Filter: 

A  special  version  of  the  original  LOWESS  filter,  dubbed  LOWESP,  was  developed  in  an  attempt  to  circum- 
vent problems  with  data  clusters  that  arise  in  many  river  and  lake  data  sets  due  to  certain  yeairs  of  intense 
sampling,  i.e.  special  surveys.  At  each  datum,  LOWESP  counts  the  number  of  observations,  say  m,  available 
within  the  specified  USI.  If  m  exceeds  the  number  of  nearest  neighbours  m  determined  from  the  choice  of 
sample  rate  and  USI,  LOWESP  performs  local  regression  with  m.  This  modification  effectively  increases  the 
smoothing  power  through  modest  data  clusters.  It  is  most  useful  when  year-to-year  smoothing  is  of  interest. 

As  per  Table  1,  LOWESP  smoothing  is  obtained  by  specifying  Filter  2. 

5.2.4  Hybrid  Filters: 

Experience  has  shown  that  there  is  no  universally  perfect  filter  capable  of  handling  all  the  potential  data  con- 
figurations presented  by  surface  water  quality  data  sets,  particularly  when  year-to-yeair  behaviour  (USI=1)  is 
of  interest.  Using  LOWESS,  LOWESP  and  Tukey's  filters,  various  hybrid  filters  (Filters  3-7  in  Table  1)  have 
been  constructed  based  on  Tukey's  notion  of  repetitive  filtering.  Filters  3  and  4  merely  polish  the  respective 
standard  LOWESS  (Filter  1)  and  LOWESP  (Filter  2)  output  with  additional  filtering  by  the  respective  par- 
ent algorithms  without  robustness  iterations  which  yield  little  benefit  when  smoothing  the  output  of  previous 
smoothing  operations.  Filters  5-7  are  somewhat  ad  hoc  constructions  that  have  in  certain  instances  been 
found  useful  for  unusual  data  configurations,  particularly  for  obtaining  a  smooth  transition  through  data 
clusters  interspersed  with  irregular  gaps  (see  conductance  example  below).  Filters  5-7  incorporate  alternate 
ways  of  smoothing  the  series  as  an  evenly  spaced  series  which  to  a  degree  circumvents  problems  created  by 
irreguleir  gaps  in  the  actual  record. 

Experience  suggests  the  following  advice  to  users: 

•  For  modestly  sized  data  series  (n  <  300),  the  hybrid  filters  can  improve  performance  through  data  series 
with  modest  data  clusters  intersperse  with  irregular  gaps  in  the  record. 

•  For  very  large  data  sets  (n  >  500),  repetitive  smoothing  can  be  excessively  time  consuming.  Use  of  the 
basic  filters  (Filter  1  and  Filter  2)  is  recommended. 

•  For  USI  >  2  there  may  be  little  advantage  to  employing  hybrid  filters. 

5.2.5  EXAMPLES:  LOWESS  and  LOWESS  derived  filters: 

The  9  robust  smoothing  filters,  including  the  Tukey  filter,  currently  available  in  programs  TRl  and  TR2  are 
listed  in  Table  1. 

5.2.5.1  Comparative  Filter  Performance:  Year-to-Year  Smoothing 

The  sequence  of  Figures  4.1-4.9  illustrate  the  behaviour  of  the  9  robust  smoothing  filters  for  a  particular 
conductance  series  where  the  USI  was  set  to  1  year,  in  order  to  emphasize  year-to-year  shifts  emd  the  median 
annual  sampling  interval  was  selected.  For  comparison.  Figure  4.10  presents  the  trend  line  developed  with 
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the  stepwise  mean  annucil  filter.  This  conductance  series  exemplifies  the  problems  encountered  with  data 
clusters  intersi)ersed  by  irregular  gaps  that  was  noted  above. 

When  examining  these  plots  consider  the  pre-  and  post-1977  behaviour  of  the  trend  lines  merits  scrutiny. 
The  rather  jagged  behaviour  of  the  pre-1977  trend  smooth  obtained  with  standard  LOWESS  is  evident  on 
Figure  4.1.  Filter  2,  i.e.  standard  LOWESP,  offers  slight  improvement  as  per  Figure  4.2.  The  repetitive 
smoothing  of  Filters  3  and  4  provide  negligible  improvement  with  the  pre-1977  record,  but  noticeably  refine 
the  trend  lines  post-1977  as  in  Figures  4.2-4.3.  Filter  5  produces  a  smooth  trend  line  over  the  entire  record 
as  per  Figure  4.5.  Filter  5  output  most  closely  matches  the  output  of  the  stepwise  filter  shown  in  Figure 
4.10  and  esthetically  produces  the  best  result  for  the  conductance  series.  Compound  Filters  6  and  7  also 
have  some  difficulty  with  the  pre-1977  records,  while  less  than  satisfactory  results  cire  obtained  with  Filters 
8  and  9  which  are  really  intended  only  for  exploratory  use  with  program  TRl. 

5.2.5.2  Comparative  Filter  Performance:  USI  =  2  years. 

The  sequence  of  Figures  5.1-5.8  illustrates  that  the  effects  of  doubling  the  unit  smoothing  interval,  i.e.  setting 
USI  =  2  years,  are  comparable  to  the  effects  of  repetitive  filtering  obtained  with  Filter  5  and  USI  =  1  year. 
When  choosing  USI  >  2,  user  should  probably  confine  filter  choices  to  Filters  1-4.  Filter  8,  the  Tukey  filter, 
always  produces  output  identical  to  Figure  4.8  as  it  is  unaffected  by  the  choice  of  sampling  interval  and  USI. 

5.2.5.3  Comparative  Performance:  Filter  4,  USI  =  3,  5,  10  &:  16.6  years. 

The  sequence  Figure  6.1-6.4,  shows  the  output  of  Filter  4  under  increasing  USI  of  3.  5,  10  and  16.6  years. 
Figures  4.4  and  5.4  should  also  be  examined  for  the  comparable  results  with  USI  =  1  and  2.  For  data  series 
like  the  conductance  example  with  considerable  year-to-year  movements,  trend  is  probably  best  represent  by 
choosing  USI  =  1  or  2.  For  series  with  minor  year-to-year  fluctuation  and  gradual  long  term  changes,  higher 
USI  may  be  in  order,  to  overcome  the  visually  obscuring  influence  of  marginal  year-to-year  fluctuations. 

The  nearly  straight  lines  indicating  downward  trend  with  USI  =  5  and  10  years,  are  graphical  analogies  of 
the  results  that  can  be  expected  of  formal  trend  tests  that  measure  net  trend  embedded  in  a  series.  The 
conductance  series  with  its  trend  reversal  since  1985  provides  a  good  example  of  why  reliance  on  formed  tests 
alone  should  be  avoided. 

5.2.5.4  Exploratory  Analysis  with  TRl. 

The  primary  uses  of  program  TRl  are: 

1.  trend  analysis  when  seasonal  effects  are  minimal  (see  next  section) 

2.  exploration  of  unfamiliar  series. 

Trend  plots  of  modestly  large  series,  e.g.  n  >  100,  look  much  like  those  produced  by  program  TR2.  Figure 
7.1  illustrates  the  usage  of  TRl  for  exploring  an  unfamiliar  series.  Choosing  USI  =  .5  year,  the  number  of 
nearest  neighbours  m  is  reduced  to  7  and  the  recurrent  oscillations  in  the  trend  line  indicate  that  significant 
seasonality  is  likely  embedded  in  the  series. 

Figure  7.2  illustrates  a  nonconventional  usage  of  TRl  that  might  be  of  use  in  some  instances.  The  plot 
shows  a  the  result  of  running  TRl  on  a  very  small  data  series  (13  observations)  which  was  obtained  during 
a  brief  period  of  5  months.  By  setting  USI  =  .1  year,  LOWESP  smoothing  is  performed  with  only  3  points, 
effectively  leading  to  a  'join  the  dots"  plot.  Far  more  efficient  software  could  be  produced  to  generate  such 
plots;  however,  as  Figure  7.2  demonstrates,  program  TRl  can  be  adopted  for  this  purpose  when  required. 

5.3  Seasonal  Smoothing  Filter: 

5.3.1  Seasonal  Time: 

Conventional  notions  of  seasons  are  generally  linked  to  arbitrarily  fixed  quarterly  or  monthly  time  periods. 
In  the  following  discussion  it  is  best  to  think  of  seasonality  in  terms  of  annually  recurrent  periodic  effects 
that  might  be  modelled  by  a  smooth  continuous  function  as  is  often  achieved  by  fitting  the  seasonal  as  the 
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sum  of  the  significant  Fourier  harmonics,  i.e.  periodic  trigonometric  functions.  To  fit  a  smooth  seasonal 
function  demands  an  appropriate  continuous  representation  of  'seasonal  time',  i.e.  time  within  the  year. 
Given  the  eccentricities  of  calendar  dates,  seasonal  time  was  handled  by  the  following  decimal  convention 
which  considerably  facilitates  plotting. 

As  data  are  read  into  the  software,  the  month,  day  and  time  of  sampling  are  converted  into  the  fraction 
of  a  complete  year  that  has  elapsed  since  January  1  0:00  hours  (equivalently  December  31  24:00  hours), 
e.g.  .5  year  =  July  2  12:00  EST  (Eastern  Standard  Time),  for  a  365  day  year.  The  date  conversion  routine 
accounts  for  leap  years.  The  chronological  time  is  then  given  by  the  yecir  plus  the  fraction,  e.g.  85.5  =  July 
2  12:00,  1985.  Generally  MOE  sample  times  are  not  recorded  in  EST  throughout  the  year,  but  rather  shift 
to  Daylight  Savings  Time  (DST)  when  appropriate.  This  introduces  a  small  discrepancy  in  seasonal  times, 
i.e.  samples  are  shifted  forward  of  EST  by  .000114  year  during  DST  periods;  however,  the  effects  on  the 
seasonal  function  produced  by  the  following  routines  is  generally  imperceptible. 

Table  2  gives  a  seasonal  date  conversion  chart  for  a  standard  365  day  year. 

5.3.2  Moving  Seasonal  Median  Filter: 

Seasonal  smoothing  is  performed  to  obtain  the  seasonal  term  S,  required  in  the  seasonal  decomposition 
routines  of  programs  TRX  and  TR2.  Seasonal  smoothing  is  accomplished  identically  in  both  programs  with 
a  moving  seasonal  median  [MSM]  filter  that  derives  from  the  data  density  plots  of  Chambers  et  al.  (1983) 
and  Tukey's  (1977)  moving  median  filters.  Eff"ectively  the  MSM  filter  amounts  to  a  'running  bin'  smoother. 
Like  the  robust  smoothers  employed  to  define  the  trend  term  T,  in  programs  TRl  and  TR2,  the  MSM  filter 
is  'quasi-continuous'  as  the  seasonal  function  is  defined  only  at  certain  seasonal  time  points  and  'joining  the 
dots'  produces  the  visual  effect  of  a  continuous  function.  Seasonal  values  required  at  other  times  are  readily 
obtained  by  linear  interpolation  to  good  approximation. 

The  MSM  seasonal  smoother  is  specifically  designed  to: 

•  counter  the  influence  of  outliers 

•  counter  the  influence  of  seasonal  sampling  bias 

•  mimic  the  appearance  of  a  smooth  continuous  se£isonal  function  such  as  might  be  derived  by  fitting 
Fourier  harmonics. 

Systematic  seasonal  sampling  bias  arises  frequently  at  sites  where  extra  sampling  is  not  evenly  distributed 
throughout  the  year,  e.g.  river  sites  where  extra  spring  runoff  samples  are  obtained  or  lake  sites  which  are 
heavily  sampled  over  the  summer  season. 

The  year  is  divided  into  100  intervals  of  3.65  days  or  effectively  101  seasonal  time  points  rt  including  the 
beginning  and  end  of  the  year,  i.e.  tq  =  rjoi.  After  reordering  the  data  series  by  seasonal  time,  a  window 
or  bin  of  width  W  is  centred  at  the  first  point  tq.  All  data  falling  within  Tq  ±  W/2  are  gathered  into  an 
array.  Given  seasonal  circularity,  for  tq  —  W/2,  data  are  those  at  the  end  of  the  year.  If  the  number  of 
observations  within  the  window  i  exceed  m,  a  parameter  specified  by  the  user,  the  median  (all  quartiles 
on  the  final  iteration)  is  determined.  When  (  <  fh,  the  window  W  is  expanded  in  small  increments  until 
sufficient  data  are  obtained  such  that  n  >  m  before  the  median  is  determined.  The  window  is  then  reset  to  W 
and  moved  forward  to  t^  and  the  process  repeated.  Essentially  the  sequence  of  101  medians  &ie  computed 
on  overlapping  data  sets  such  that  only  a  small  fraction  of  the  data  change  at  each  r^.  Two  rounds  of 
basic  LOWESS  smoothing  with  10  nearest  neighbours  and  0  robustness  iterations  are  applied  to  polish  the 
sequence  of  101  medians,  once  again  taJcing  advantage  of  the  circularity  of  the  seasoned  data  series.  Finally 
the  100  distinct  smoothed  medians  are  summed  and  adjusted  to  yield  a  sum  of  0.  Experience  shows  that  a 
slight  adjustment  is  needed  when  sesisonally  biased  sampling  records  axe  analyzed  to  assure  that  the  resulting 
relative  seasonal  function  is  perfectly  centred  about  0. 

5.3.3  Seasonal  Smoothing  Control  Parameters: 

The  two  seasonal  smoothing  control  parameters  aje: 
1.  m  =  required  minimum  number  of  data  falling  within  a  seasonal  window  before  computation  of  statistics 
is  performed 
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2.  h  =  number  of  seasons  per  year  such  that  seasonal  smoothing  window  width  W  =  l/h 

The  goal  in  selecting  in  and  h  is  to  achieve  as  sharp  a  definition  of  a  smoothly  varying  seasonal  variation 
as  data  permit  as  economically  as  possible.  The  selection  procedure  remains  a  somewhat  ad  hoc,  heuristic 
process.  Generally  it  is  desirable  to  determine  values  that  are  appropriate  throughout  most  of  the  yecir. 
Practice  suggests  that  when  the  seasonal  variabihty  embedded  in  the  data  series  is  reasonably  strong,  a 
reasonable  degree  of  latitude  is  possible  in  the  selection  of  fh  and  h  with  negligible  effect  on  the  results. 

The  primary  purpose  of  parameter  m  is  to  serve  as  a  fail-safe  mechanism  that  triggers  seaisonal  window 
expansion  to  take  over  through  any  thin  spots  in  the  seasonal  sampling  distribution,  e.g.  at  most  PWQMN 
river  sites,  December-January  is  typically  a  poorly  sampled  time  of  the  year.  The  mechanism  assures  a 
smooth  transition  of  the  resulting  seasonal  function  through  poorly  sampled  times  of  the  year  during  which 
the  resulting  seasonal  function  is  essentially  a  projection  of  the  seasonal  behaviour  exhibited  to  the  well 
sampled  periods  to  the  right  and  left  of  the  poorly  sampled  period. 

The  purpose  of  parameter  h  is  to  determine  the  size  of  the  seasonal  smoothing  window  W  =  l/h.  It  is 
simplest  to  envision  h  as  the  number  of  feasible  distinct  seasons  that  might  be  obtained  if  the  n  of  the  series 
were  spread  evenly  throughout  the  year.  Practice  has  shown  that  if  there  are  sufficient  data,  sharp  seasonal 
definition  is  obtained  by  choosing  h  =  12. 

As  the  process  involves  some  trial  and  error,  experiment  with  alternative  values  of  fh  and  h  using  progrcmi 
TRX  for  quick  results,  until  acceptable  values  are  found.  If  the  suggested  choices  for  fh  are  employed,  the 
associated  fail-safe  mechanism  and  the  two  rounds  of  LOWESS  smoothing  generally  maintain  the  resulting 
seasonal  function  within  acceptable  bounds.  Diagnostic  programs  to  obtain  preliminary  cissessments  of 
seasonal  sampling  bias  and  guide  the  selection  of  fh  and  h  are  under  preparation.  For  the  present  the 
following  rules  are  suggested: 

•  Keeping  in  mind  that  the  median  grows  less  reliable  as  sample  size  decreases,  pajcimeter  m,  in  general, 
should  never  be  set  <  10.  With  modest  data  series,  i.e.  n  =  150-250,  and  some  seeisoucil  sampling  bijis 
present,  values  of  m  =  11-13  were  found  satisfactory  for  many  of  the  examples  shown  in  this  guide. 

•  Begin  with  a  value  of  /i  =  12  (corresponding  to  365/12  s:  30  day  seasons),  as  this  usually  yields  sharp 
seasonal  definition  if  sufficient  data  are  available. 

•  Assuming  data  are  evenly  distributed  across  the  year,  estimate  the  expected  number  of  observations  as 

e  =  n/h. 

•  U  £  <  fh,  decrease  the  number  of  seasons  h  (equivalently  increasing  window  size  W  =  l//i)  until  a 
suitable  value  of  (  is  obtained.  Generally,  integral  values  in  the  range  /»  =  4, 12  are  sufficient  for  most 
series  of  modest  size.  Even  with  fairly  small  data  sets,  e.g.  n  =  90,  perceptible  seasoned  variation  may 
be  evident  at  /)  =  4  ss  91  day  or  3  month  sccisons  or  /i  =  6  «  61  day  or  2  month  seasons. 

•  For  large  data  series,  e.g.  n  >  300,  users  may  consider  /i  =  17  ss  21  day  or  3  week  seasons  or  /i  =  26  ss  14 
day  or  2  week  seasons.  When  working  with  large  data  sets,  if  ^  =  n/h  >  fh,  users  should  reset  m  ss  £. 

5.3.4  EXAMPLES:  Seasonal  Smoothing: 

The  following  series  of  examples  were  culled  from  test  data  sets  and  are  selected  to  illustrate  various  types 
of  seasonal  trends  and  other  points  regarding  seasonal  functions  developed  by  the  TRENDS  software. 

5.3.4.1  Clearwater  Lake  Examples: 

Clearwater  Lake  is  a  small  Shield  lake  located  in  the  Sudbury  area.  The  chemistry  and  seasonal  dynamics 
indicated  by  the  following  examples  are  indicative  of  acidified,  dilute  waters. 

Figure  8.1        Seasonal  trends  in  DIC,  Clearwater  L. 

•  strong  seasonsility  is  suggested  concordant  with  the  exploratory  results  obtained  with  program  TRl, 
Figure  7.1. 

•  a  case  where  very  strong  seasonality  is  evident  in  a  relatively  small  (n  =  91)  data  set. 
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Figure  8.2       Seasonal  trends  in  conductance,  Clearwater  L. 

•  An  example  of  reasonably  strong  seasonality  encountered  in  the  conductance  series  that  was  employed 
in  the  Figure  4.1-4.10  sequence.  Figure  8.2  was  produced  specifically  as  the  companion  for  Figure  4.5; 
however,  Filters  1-7  all  produce  results  that  are  very  close  to  Figure  8.2 

Figure  8.3       Seasonal  trends  in  manganese,  Clearwater  L. 

•  An  example  of  modest  seasonality  that  is  unlikely  to  have  adverse  effects  on  trend  analysis  but  may 
merit  notice  from  a  scientific  perspective. 

•  An  example  of  scaling  problems  with  PLOT88.  Notice  the  upper  40%  of  the  graph  sp«ice  is  blank, 
flattening  the  appearance  of  the  seasonal  function  and  artificially  creating  the  impression  that  the 
seasonal  my  be  less  significant  than  it  is. 

Figure  8.4       Seasonal  trends  in  potassium,  Clearwater  L. 

•  An  example  of  a  weak  secisonal  function  that  has  little  or  no  eff'ect  on  the  results  of  trend  analysis 
but  may  yet  be  significant,  i.e.  nonrandom.  When  percentage  variation  contributed  by  the  seasonal  is 
<10%  after  outlier  deletion  (see  PART  II:  5.2),  the  seasonal  function  estimated  by  the  software  is  likely 
insignificant. 

Figure  8.5a— d       Seasonal  trends  in  aluminum,  Clearwater  L: 
Vtu-ying  seasonal  smoothing  control  parameters. 

•  Figures  8.5a-8.5c  show  the  effects  of  varying  seasonal  smoothing  control  parameters  fh  and  /i  on  an 
aluminum  data  set  of  modest  size.  The  choices  were: 

8.5a  m  =  11  and  h  =  10,  i.e.  W  =  36.5  days 
8.5c  fn  =  13  and  /j  =  8,  i.e.  W  =  45.6  days 
8.5d  m  =  15  and  /i  =  6,  i.e.  W  =  60.8  days 

•  In  moving  from  Figure  8.5a  to  8.5c,  the  estimated  seasonal  function  does  become  smoother  as  more 
data  are  considered  at  each  of  the  100  points  at  which  the  seasonal  is  estimated. 

•  The  seasonal  data  series  is  characterized  by  considerable  scatter,  but  all  results  suggest  that  there  may 
be  net  seasonality  of  minor  significance  embedded  in  the  data  series.  There  is  probably  no  best  choice. 
It  is  best  not  to  focus  on  specific  features,  but  rather  to  note  general  shape  and  degree  of  scatter. 

•  Figure  8.5d,  the  trend  plot  produced  with  the  choices  for  Figure  8.5a,  shows  an  dominant  downward 
trend.  All  choices  in  fcict  yield  the  same  trend  result  as  would  an  analysis  with  program  TRl. 

Figure  8.6a-d       Seasonal  trends  in  chloride,  Clearwater  L: 

Contrast  of  seasonal  estimated  via  TRX  against  TR2  under  constant  trend. 

•  Figures  8.6a-8.6b  show  respectively  the  trend  and  seasonal  plots  produced  for  the  Clearwater  L  chloride 
series  with  program  TRX. 

•  Figures  8.6c-8.6d  show  respectively  the  trend  and  seasonal  plots  produced  for  the  Clearwater  L  chloride 
series  with  program  TR2. 

•  Steady  upward  trend  is  the  dominant  feature  of  the  trend  plots  produced  by  either  programs;  however, 
the  seasonal  plots  show  somewhat  different  characteristics. 

•  Under  a  steady  upward  or  or  downward  trend  the  seasonal  plots  produced  by  program 
TRX  will  distort  the  shape  and  other  features  of  the  seasoned  function.  When  estimating 
trend  as  a  constant  level  throughout  each  year,  in  the  case  of  the  chloride  data  with  steady  upward 
trend,  the  level  of  the  seasonal  during  the  early  pcirt  of  the  year  is  too  low  eind  during  the  latter  part  of 
the  year  is  too  high. 

•  To  investigate  the  form  of  seasonal  characteristics  under  steady  unidirectional  trend  it  is 
desirable  to  employ  program  TR2. 

•  When  there  is  either  little  trend  or  year-to-year  upward  trend  more  or  less  equals  year-to-year  downward 
trend  over  the  time  horizon  of  a  series,  there  should  be  little  difference  observed  between  the  seasonal 
chairacteristics  estimated  by  TRX  and  TR2. 
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5.3.4.2  Moira  River  Arsenic  Examples: 

Varying  seasonal  smoothing  control  parameters. 

The  upper  reaches  of  the  Moira  River  are  contaminanted  with  Arsenic  as  the  result  of  long  abandonned 
gold  mining  operations.  As  the  following  sequence  of  figures  illustrate  that  with  a  large  data  set,  i.e.  400 
observations,  and  strong  seasonality,  considerable  latitude  is  possible  in  the  selection  of  control  parameters 
m  and  h  with  only  minimal  effects  on  the  resulting  seasonal.  The  choices  were: 

(9.a)  m  =  15  and  /i  =  26,  i.e.  W  %  14  days 

(9.b)  m  =  35  and  /i  =  12,  i.e.  W  k  30  days 

(9.0)  m  =  70  and  /i  =  6,  i.e.  W  k  61  days 

Figures  9.a-9.c  speak  for  themselves.  Figure  9.d  shows  the  historical  trend  plot  generated  as  a  companion 
to  Figure  9.b.  The  choice  of  smoothing  filter  control  parameters  had  no  perceptible  influence  on  trend  plots. 
Figure  9.d  is  another  example  of  the  scaling  problems  occasionally  encountered  with  PLOT88  software. 

5.3.4.3  Nitrate  Seasonality  Examples: 

Water  quality  variables  nitrates,  i.e.  NOj,  and  total  nitrates,  i.e.  NOj  +  NOj,  frequentlv  show  strong 
seasonal  behaviour  as  in  the  following  examples. 

Figure  10.1        Seasonal  trends  in  total  nitrates,  Trent  R. 

•  extremely  strong  unimodal  seasonality  is  evident  in  a  very  large  data  set. 

•  the  Trent  R  represents  a  mixture  of  Shield  and  Lowland  drainage. 

Figure  10.2        Seasonal  trends  in  nitrates.  Plastic  L. 

•  Plastic  L  is  a  dilute,  acidified  lake  in  central  Ontario  Shield  area. 

•  strong  bimodal  seasonality  is  evident  in  the  Plastic  L  inflows. 
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PART  II:  RUNNING  TRENDS  PROGRAMS: 


0.0  HARDWARE  REQUIREMENTS: 

The  current  edition  of  the  TRENDS  programs  is  designed  for  implementation  on  a  386  IBM 
PC  or  PC— clone  with  a  80387  math  coprocessor  operating  under  DOS  3.31  or  higher  and 
connected  to  an  HP  LaserJet  Series  II  laser  printer. 

The  laser  printer  should  have  at  least  4  megabytes  of  extended  memory  in  order  to  hold  the 
graphics  images  produced  by  TRENDS  programs. 

1.0  PROGRAM  EXECUTION: 

Programs  are  activated  by  typing  the  program  name  (without  suffix),  e.g.  to  run  TRX,  simply  type  TRX 
and  hit  the  ENTER  key.  On  entry  of  the  desired  program  name,  the  computer  responds: 

TYPE  CONTROL  FILE  NAME: 
whereupon  the  user  must  enter  the  name  of  the  file  containing  the  program  control  information  for  that 
particular  run,  e.g.  GRAND. CNT.  On  hitting  the  enter  key,  the  program  proceeds  to  execute  to  completion 
unless: 

1.  errors  with  the  control  information  or  the  input  data  set  are  encountered 

2.  Itiser  printer  problems  are  encountered. 

During  the  development  of  the  current  edition  of  these  programs,  numerous  checks  for  illegal  control  param- 
eters and  input  data  indicators  were  incorpwrated  into  the  software.  In  some  cases  the  programs  will  stop. 
Whenever  a  program  crsishes,  users  should  immediately  check  the  report  file  for  error  messages.  Users  may 
encounter  program  stoppages  due  to  difficulties  as  yet  unencountered  by  the  programs'  author. 

2.0  PROGRAM  SIZE  LIMITS: 

Due  to  variations  in  program  capabilities,  the  size  limitations  on  the  input  data  files  handled  by  the  respective 
programs  vary  as  indicated  below.  The  limits  given  below  cire  for  the  most  compact  and  fastest  versions  of 
the  programs.  These  programs  should  handle  90%  of  water  quality  data  series  encountered  in  practice. 

TRl:  maximum  series  length  =  1000  elements 

TRX:  maximum  series  length  =  900  elements 

TR2:  mciximum  series  length  =  800  elements 

Very  Large  Data  Sets 

Special  huge  versions  of  the  programs  have  been  compiled  to  handle  very  large  water  quality  time  series. 
The  huge  versions  of  these  programs  were  achieved  by  sacrificing  error  diagnostics,  some  statistical  routines 
cind  slower  execution  times.  All  data  input  and  program  control  requirements  are  the  same  as  for  the  TR* 
series.  These  programs  have  not  been  extensively  tested. 

THl:  maximum  series  length  =  3000  elements 
THX:  maximum  series  length  =  2500  elements 
TH2:  maximum  series  length  =  2000  elements 

3.0  INPUT  DATA  REQUIREMENTS: 

The  input  data  file  formats  are  identical  for  all  3  programs.  Input  data  must  be  in  the  form  of  an  ASCII  file 
formatted  according  to  standard  FORTRAN  conventions  a.s  per  the  example  given  in  Table  3. 
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3.1  Line  1:   Parameter  i:  Site  Identifiers 

The  first  line  of  the  input  data  file  should  contain: 

•  cin  abbreviated  water  quality  parameter  name  in  the  first  6  columns,  e.g.  'As'  (Arsenic)  in  the  Table  3 
example  or  'COND'  (conductance)  as  in  the  Table  4  example. 

•  a  12  chariicter  site  identification  code  in  columns  7-18,  e.g.  '170026013  '  in  the  Table  3  example  or 
'Clearwater  L'  as  in  the  Table  4  example. 

If  either  the  parameter  abbreviation  is  less  that  6  characters  or  the  site  identifier  is  less  than  12  characters, 
these  should  be  left  justified  to  the  first  allotted  column,  i.e.  columns  1  and  7  respectively. 

N.B.  The  parameter  name  and  site  identifier  in  line  1  of  the  input  data  file  appear  on  the  graphics.  This 
feature  was  introduced  after  observing  that  staff  hastily  generating  graphical  results  only  with  test  versions 
of  the  programs  were  forgetting  to  change  plot  headers  and  y-axis  titles  and  consequently  confusing  sites 
and  variables.  Although  not  entirely  fool  proof,  this  feature  provides  a  checking  mechanism  in  the  event  that 
program  users  change  input  data  files  without  revising  plot  header  and  y-axis  information. 

3.2  Lines  2-end  of  file:  Input  Data  Series 

The  TRENDS  programs'  input  routines  seek  the  following  three  variables: 

1.  sample  date  read  as  a  6  digit  integer  indicating  year/month/day.  e.g.  641101  for  1964  November  1, 

2.  sample  time  read  as  a  4  digit  integer  indicating  hour/minutes,  e.g.  1430  for  2:30  p.m. 

3.  sample  concentration  read  in  FORTRAN  real  fbced  format,  i.e.  Tx.y  where  x  is  the  total  width  of  the 
sample  concentration  field  and  y  <  i  is  the  number  of  digits  to  the  right  of  the  decimal. 

In  the  example  data  file  of  Table  3,  the  first  9  columns  of  lines  >  2  happens  to  be  occupied  by  the  site 
identification  code.  The  sample  date  appears  in  columns  11-16,  the  sample  time  in  columns  18-21  and  the 
arsenic  concentrations  (in  /Jg/L)  appear  in  columns  22-28.  Either  of  the  following  FORTRAN  formats  would 
lead  to  data  being  correctly  read  by  the  trend  programs: 

•  (10X,I6,I5,F7.0) 

•  (10X,I6,1X,I4,F7.0) 

where  the  'lOX'  skips  the  reader  over  the  1st  10  columns  of  the  line  and  'I;'  indicates  an  integer  field  of 
length  z. 

Sample  times  are  of  minor  importance  with  sparse,  e.g.  monthly,  sampling;  however,  when  records  contadn 
instances  of  more  than  1  sample  per  day,  sample  times  are  important  for  distinguishing  between  samples 
obtained  on  the  same  day.  In  SIS,  missing  sample  times  are  given  the  vaJue  zero.  The  trend  softwcire  converts 
zero  times  to  12:00  noon  since  virtually  no  PWQMN  samples  have  ever  been  collected  at  midnight. 

N.B.  If  data  records  exist  for  which  there  are  no  sample  times  available,  the  input  file  should  be  prepared 
with  at  least  1  column  of  blanks  separating  the  sample  date  and  the  sample  concentration.  If  there  were 
1  blank  column  between  date  and  concentration,  the  input  format  could  be  given  as  (I6,Il,Fx.y)  by  which 
sample  times  would  be  read  as  an  II  column  of  blanks  that  are  interpreted  as  zeros  and  converted  to  12:00 
noon. 

For  convenience  of  numerical  manipulation  and  plotting,  the  software  internally  determines  a  real  sample 
decimal  date,  e.g.  YR.xxxxx  where  the  YR  is  the  year,  .e.g.  75,  and  .xxxxx  is  the  hours,  minutes,  day 
jind  month  of  the  sample  expressed  as  a  decimal  fraction  of  the  year.  Leap  years  are  identified  cind  haindled 
correctly. 

N.B.  When  entering  input  formats  in  the  control  (.CNT)  files,  the  parentheses  must  be  included. 

N.B.  The  trend  programs  were  designed  to  work  with  concentration  data  having  no  more  than  3  digits  to 
the  right  of  the  decimal.  The  remark  applies  specifically  to  cases  like  Cadmium  which  if  retrieved  with  SIS 
default  options  gives  typical  sample  concentrations  as  .OOOi  mg/L  or  .OOxy  mg/L.  Data  of  this  type  should 
always  be  converted  to  more  appropriate  units  as  pg/L  or  ng/L.  For  PWQMN  river  data,  this  can  always 
be  done  by  specifying  appropriate  unit  conversions  when  retrieving  data  from  SIS.  If  a  data  file  has  been 
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prepared  with  concentration  data  have  >  4  digits  to  the  right  of  the  decimal,  a  unit  conversion  feature  exists 
in  the  .CNT  files  that  would  allow  such  input  data  to  be  converted  to  units  with  <  3  digits  to  the  right  of 
the  decimal  (see  PART  II,  4.2  Line  4).  Note  that  for  data  like  those  presented  in  Table  4,  the  trailing  digits 
to  the  right  of  the  decimal  have  no  effect  on  the  analysis. 

N.B.  Current  versions  of  the  trend  programs  cannot  read  the  ubiquitous  alphcinumeric  remark  codes  that 
accompany  numerical  values  of  data  retrieved  from  SIS.  The  '<'  symbol  appearing  to  the  right  of  the 
concentration  value  of  the  sample  obtained  on  date  671002  of  Table  3  is  an  example.  If  input  data  format 
is  incorrectly  specified  and  alphanumerics  are  encountered,  the  program  will  crash  and  return 
the  FORTRAN  error  message: 

•  run-time  error  F6103:  RE AD(input filename)  -  invalid  real 

3.3  Data  Matrices: 

The  input  format  requirements  make  it  possible  to  prepare  a  multi-variable  input  file  for  a  single  site  as 
shown  in  the  example  of  Table  4  where  the  date  occupies  columns  2-7,  the  time  columns  9-12  and  the  5 
variables  were  formatted  as  5F12.4.  The  5  variables  happen  to  be: 

1.  conductance  {/iS/cm  @  25C) 

2.  hydrogen  ion  concentration  (/jeq/L) 

3.  calcium  concentration  (/ieq/L) 

4.  magnesium  concentration  (^eq/L) 

5.  sodium  concentration  (/ieq/L) 

N.B.  The  missing  value  code  -999.0000  (see  PART  II:  4.2  Line  6)  has  been  entered  in  the  calcium  and 
magnesium  columns  as  data  records  for  these  variables  do  not  begin  until  1977. 

N.B.  Line  1  of  a  multi-variable  file  is  identical  to  that  of  a  single  variable  file.  When  working  with  an  input 
data  matrix,  it  is  more  convenient  to  enter  some  generic  variable  indicator  in  columns  1-6  of  the  first  line  as 
in  Table  4,  rather  than  repeatedly  enter  and  edit  the  variable  name,  rather  then 


4.0  PROGRAM  CONTROL  FILES: 

Programs  execution  features  are  controlled  by  .CNT  (or  einy  suffix  the  user  may  prefer)  control  files.  For 
consistency  and  ease  of  use,  control  files  for  the  three  main  trend  analysis  programs  are  structured  similarly, 
differing  only  by  one  or  two  lines  according  to  the  particular  requirements  of  the  program.  Sample  .CNT 
files  for  the  three  programs  TRl,  TR2  and  TRX  are  given  in  Table  5.  The  common  elements  of  the  three 
programs  are  listed  below: 

4.1  Common  Lines  of  TRl,  TR2  &  TRX  Program  Control  Files: 

1.  lines  1-7:  all  three  programs 

2.  line  8:  TRl  [1st  2  elements],  TR2  i:  TRX 

3.  line  9:  TRl  k.  TR2 

4.  TRX  (line  9)  i:  TR2  (line  10) 

5.  TRl  k.  TRX  (lines  10-17);  TR2  (lines  11-18) 

4.2  Control  File  Input  Requirements: 

Line  1:  all  3  programs;  Report  File  Name 

•  48  character  report  file  name,  left  justified  to  column  1 
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•  the  report  file  is  read  1st  because  it  is  opened  immediately  to  record  subsequent  information  in  the 
.CNT  file 

Line  2:  all  3  programs;  Input  Data  File  Name 

•  48  character  input  data  file  name,  left  justified  to  column  1 

Line  3:  all  3  programs;  Water  Quality  Variable  Name 

•  6  character  water  quality  parameter  name,  left  justified  to  column  1 

Line  4:  all  3  programs;  Input  Data  Format 

•  40  character  (including  parentheses)  input  data  format,  left  justified  to  column  1 

Line  5:  all  3  programs;  Data  Bracketing  Dates 

•  data  series  bracket  dates,  i.e.  only  samples  occurring  within  brackets  are  read  as  input 

-  lower  bracket  date  as  6  digit  integer  (year/month/day)  in  columns  2-7 

-  upper  bracket  date  as  6  digit  integer  (year/month/day)  in  columns  9-14 

Line  6:  all  3  programs;  Zero  Option  &:  Missing  Value  Codes 

•  1  digit  zero  option  code  [ZOC]  (column  1) 

-  set  to  0  to  leave  O's  as  is 

-  set  to  1  to  replace  O's  with  1/2  smallest  positive  value  in  record 

-  set  to  2  to  replace  O's  with  smallest  positive  value  in  record 

-  set  to  3  to  delete  O's 

•  missing  value  code  [MVC]  (columns  2-10);  format  F9.3 

-  N.B.  MVC  must  be  a  negative  value 

-  for  data  sets  with  legitimate  negative  values,  e.g.  Gran  alkalinity  k.  residual  series,  set  MVC  to  a 
value  smaller  than  the  smallest  member  of  series 

-  N.B.  F9.3  format  can  be  over  ridden  if  necessary,  e.g.  F9.5,  and  will  be  read  properly.  The  desired 
MVC  must  appear  in  columns  2-10  and  include  the  decimal. 

Line  7:  all  3  programs;  Multiplier/Divisor  Codes 

•  1  digit  multiplier/divisor  code  [MDC]  (column  1);  format  Al 

-  set  to  1  for  no  multiplier/divisor,  i.e.  to  leave  data  as  is 

-  set  to  M  to  multiply  input  data 

-  set  to  D  to  divide  input  data 

-  N.B.  Codes  M  i:  D  must  be  upper  case 

•  multipher/divisor  (columns  2-10);  format  F9.3 

-  if  MDC  =  1,  can  be  ignored 

-  if  MDC  =  M  or  D,  specify  numerical  value 

-  N.B.  F9.3  format  can  be  over  ridden  if  necessary,  e.g.  F9.5,  and  will  be  read  properly.  The  desired 
multipher/divisor  must  appear  in  columns  2-10  and  include  the  decimal. 

Line  8:  all  3  programs;  Log,  Run  &:  Iteration  Parameters 

•  1  digit  logio  transform  code  (column  3) 

-  set  to  0  for  no  logio  transform 

-  set  to  1  for  logio  transform 

•  #  of  program  runs  (see  PART  11:  5.2),  3  digit  integer  or  less  (preferrably),  right  justified  to  column  6 

-  set  to  2  or  3  for  trial  runs 

•  #  of  iterations  in  seasonal  decomposition  loop,  3  digit  integer  or  less  (preferrably),  right  justified  to 
column  9 

-  3  iterations  are  sufficient  for  most  work 
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Line  9:  TRl  &  TR2;  Smoothing  Filter  Parameters 

•  1  digit  smoothing  filter  identifier  (column  2) 

-  choose  from  filters  1  to  9  (PART  I:  Section  5) 

•  1  digit  sampling  interval  code  (column  4) 

-  choose  1  for  median  sampling  interval  (most  data  sets) 

-  choose  2  for  mean  sampling  interval  (seasonally  sampled  lake  k.  river  records) 

•  unit  smoothing  interval;  format  F6.2  (columns  5-10) 

-  typically  1  or  2  years  with  TR2 

Line  9  (TRX);  Line  10  (TR2);  Seasonal  Smoothing  Parameters 

•  minimum  number  of  observations  required  to  compute  statistics  within  seasonal  smoothing  windows 

-  set  to  11  or  higher  (see  PART  I:  Section  5.3  L  PART  II  Section  5.4) 

•  minimum  width  of  seasonal  smoothing  window 

-  set  as  a  fraction  of  1  year,  e.g.  12  indicates  1/12  year  (see  PART  I:  Section  5.3  k.  PART  II  Section 

5.4) 

Line  10  (TRX  &:  TRl);  Line  11  (TR2);  Plot  Codes 

•  plot  code:  1  digit  character  in  column  2 

-  set  to  upper  case  'P'  for  laser  printer  plots 

-  anything  else  fails  to  produce  laser  printer  plots 

•  plot  density  codes:  1  digit  character  in  column  4 

-  set  to  upper  case  'S'  for  single  density  plots  (75  dots  per  inch;  suitable  for  fast,  rough  work) 

-  set  to  upper  case  'D'  for  double  density  plots  (150  dots/inch;  suitable  for  most  reports) 

-  set  to  upper  case  'Q'  for  quad  density  plots  (300  dots/inch;  suitable  for  high  quality  reports  k. 
journal  articles) 

•  plot  type  codes:  1  digit  character  in  column  6 

-  set  to  '0'  for  ordinary  trend  plots  of  either  untransformed  or  log  transformed  variables 

-  set  to  '1'  for  logarithmic  scaled  plot  (only  permissible  when  log  transform  code  =  1 

•  plot  annotation  code:  1  digit  character  in  column  8 

-  set  to  upper  case  'S"  for  plot  without  right  side  annotation 

-  set  to  upper  case  'L'  for  plot  with  right  side  annotation 

Line  11  (TRX  &  TRl);  Line  12  (TR2):  Plot  Header  1 

•  up  to  56  characters,  left  justified  to  column  1 

-  appears  below  main  plot  header 

Line  12  (TRX  &:  TRl);  Line  13  (TR2):  Plot  Header  2 

•  up  to  56  characters,  left  justified  to  column  1 

-  appears  below  plot  header  1 

Line  13  (TRX  &  TRl);  Line  14  (TR2):  Current  Date 

•  up  to  20  characters,  left  justified  to  column  1 

-  e.g.  May  31  1991 

-  appears  in  upper  right  corner 

Line  14  (TRX  L:  TRl);  Line  15  (TR2):  I'-axis  Title 

•  up  to  30  characters,  left  justified  to  column  1 

-  if  log  transform  is  specified,  program  automatically  supplies  'Logio'  prefix  where  necessary 

Line  15  (TRX  &  TRl);  Line  16  (TR2):  Report  Code 
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•  2  digit  character  code  left  justified  to  column  1 

•  specify  'SR'  for  short  report 

-  short  reports  do  not  include  input  data  list  and  skip  correlation  calculations  on  the  final  results 

-  short  reports  are  recommended 

•  specify  'LR'  for  long  report 

-  long  reports  include  input  data  list  and  correlation  calculations  on  the  final  results 

-  the  correlations  are  experimental  and  may  be  extremely  time  consuming 

-  long  reports  are  not  recommended 

Line  16  (TRX  &:  TRl);  Line  17  (TR2):  Output  File  Code 

•  1  digit  code  left  justified  to  column  1 

-  set  to  '0'  for  no  output  files 

-  set  to  any  non-'O'  1  digit  character  for  output  files 

Line  17  (TRX  &:  TRl);  Line  18  (TR2):  Trend  Output  File  Name 

•  see  Section  6.6 

Line  18  (TRX);  Line  19  (TR2):  Seasonal  Output  File  Name 

•  see  Section  6.6 


5.0  TRENDS  Primary  Program  Control  Parameters: 

There  are  3  fundamental  program  control  options  that  govern  the  execution  of  TRENDS  programs: 

1.  log  code 

2.  number  of  runs 

3.  number  of  iterations  (not  necessary  for  program  TRl). 

These  are  entered  in  Line  8  of  the  program  control  files  of  all  three  programs.  These  are  discussed  separately 
because  of  their  primary  importance.  The  remaining  control  parameters  are  discussed  in  PART  II:  Section 
6. 

5.1  Data  Transformation: 

The  present  TRENDS  software  permits  only  two  transforms: 

•  Z,=  X,  (null) 

•  Zi  =  logio  A', 

More  precise  transformations  other  than  logarithmic  from  the  broad  Box-Cox  family  of  power  trcinsforms 
might  be  necessary  in  certain  cases  when  more  formal  modelling  exercises  are  undertaken.  Because  the 
statistical  methods  employed  in  TRENDS  programs  are  highly  robust,  even  if  neither  the  null  nor  the  logjg 
transforms  are  the  most  appropriate,  at  least  one  of  one  of  the  two  available  choices  should  lead  to  satisfactory 
graphical  results. 

Many  variables,  in  particular  those  exhibiting  positive  correlation  with  flow,  are  positively  skewed,  i.e.  many 
observations  cire  close  to  the  lower  bound  and  relatively  few  lie  one  or  more  orders  of  magnitude  above  the 
lower  values.  For  positively  skewed  data,  logarithmic  transforms  commonly  lead  to  data  distributions  that 
more  closely  approximate  the  normal.  Even  if  raw  data  are  not  necessarily  positively  skewed,  when  data 
rcinge  over  more  than  1  order  of  magnitude,  a  logarithmic  transformation  may  be  advisable  as  means  of 
obtaining  good  graphical  results.  In  time  series  analysis,  it  is  also  desirable  that  local  variabihty  exhibited 
by  the  data  series  be  relatively  constant  over  the  time  horizon  of  the  series.  Again  a  log  transform  may  help 
achieve  this  objective.  However,  this  rule  may  not  necessarily  be  apropos  if  the  wide  variation  in  the  data 
series  is  attributable  to: 

•  only  a  few  outliers 
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•  significant  trend  in  the  presence  of  limited  relatively  constant  local  variation. 

For  the  case  where  wide  variability  is  caused  by  only  a  few  outliers,  it  is  preferable  to  leave  data  as  is  and  rely 
on  the  runs  feature  described  below  to  eliminate  the  troublesome  extreme  values.  Marginal  cases  do  arise; 
however,  exjjerience  suggests  that  the  results  of  either  a  null  or  a  logio  transform  should  be  comparable. 

When  the  logjo  transform  is  invoked,  it  should  be  kept  in  mind  that  models  (1-3)  aie  linear  additive  on 
the  logarithms  and  that  care  must  be  taken  with  the  antilog  transformation  and  the  interpretation  of  the 
results.  With  the  seasonal  axijustment  model  (3),  the  relative  seeisonal  comp>onent  S,  must  be  added  to  some 
trend  level,  e.g.  the  time-weighted  mean,  before  antilogs  are  taken,  in  order  to  have  meaning  on  the  natural 
scjde.  The  seasonal  amplitude  is  constant  only  on  the  logarithms.  On  the  natural  scale  seasonal  periodicity 
varies  according  to  whatever  level  T  on  which  the  relative  seasonal  function  has  been  superimp>osed.  For 
the  written  report  file,  when  data  have  been  logjo  transformed,  the  programs  TRX  and  TR2  determine  the 
seasonal  amplitude  at  three  trend  levels: 

1.  minimum  annual  level 

2.  overall  geometric  mean  trend  level 

3.  maximum  annual  level 

For  program  TRX,  the  minimum  and  maximum  annual  levels  are  taken  as  the  years  with  the  smallest  and 
largest  medians  as  determined  on  the  de-seasonalized  data,  while  in  program  TR2,  the  respective  levels 
are  determined  as  the  yearly  time-weighted  means  of  trend  estimates  T,  occurring  within  the  given  years. 
For  TRX,  the  overall  geometric  mean  trend  is  the  mean  of  previously  cited  yearly  medians.  For  TR2,  the 
overall  geometric  mean  level  is  determined  from  the  grand  time-weighted  mean  trend,  calculated  with  trend 
estimates  T,  as  determined  for  the  logio  transformed  series.  Finally,  the  results  are  antilogged  to  permit 
users  to  examine  the  range  of  seasonal  variation  at  the  three  trend  levels  back  in  the  natural  scale. 

5.2  Runs:  Eliminating  Extreme  Values: 

As  the  final  step  of  a  complete  analysis  or  "run'  for  a  given  data  set,  the  three  programs  TRX,  TRl  L  TR2 
perform  an  analysis  of  the  model  residuals,  e.g. 

£,  =  X,-T„         (TRl) 

e,  =  A'.  -  T,  -  Si,         (TRX,  TR2) 

in  order  to  identify  extreme  outliers  which  are  clcissified  according  to  various  outlier  identification  criteria. 
At  the  final  step  of  the  residual  analysis,  extreme  data  identified  as  'Far  Outliers'  or  'Far  Out  Points* 
according  to  Tukey's  (1977)  criterion  are  deleted  from  the  data  set.  For  a  normal  distribution,  Tukey's 
'Outer  Fences"  coincide  with  Qso  ±  4.72cr,  where  Qso  is  the  median  and  a  is  the  standard  deviation.  Of 
several  available  choices,  this  outlier  identification  criterion  is  the  most  conservative  in  the  sense  that  it  is 
the  widest  and  leads  to  the  deletion  of  the  fewest  extrema.  Indeed  Tukey  Outer  Fences  bound  the  central 
99.999767o  of  the  normal  distribution. 

The  TRENDS  programs  maintain  a  running  count  of  completed  'runs'  and  after  a  run  is  complete,  the  count 
is  checked  against  a  user  specified  run  limit  parameter  in  the  main  program  control  line  of  the  .CNT  file.  If 
the  completed  runs  are  less  than  the  user  specified  run  limit,  the  analysis  is  repeated  on  the  reduced  data 
set  free  of  the  Far  Outliers  identified  at  the  end  of  the  previous  run.  If  no  far  outliers  were  identified  on  the 
previous  run,  the  program  stops. 

In  general  there  are  no  hard  rules  for  specifying  the  number  of  runs.  A  minimum  of  3  is  recommended; 
however,  it  is  the  author's  preference  to  specify  a  high  limit,  e.g.  10,  and  let  the  programs  run  to  completion. 
For  reporting  purposes,  analyses  may  be  repeated  to  some  satisfactory  level,  e.g  2  or  3.  Frequently  the  first 
round  of  outlier  deletion  is  essential  to  eliminate  extreme  distortion  of  the  graphics  y-axis  scales.  Users 
should  not  be  unduly  concerned  with  the  number  of  Far  Outliers  identified  if  these  constitute  <  5%  of  the 
members  of  the  original  data  set.  If  analyses  are  being  performed  on  untransformed  data  and  the  number 
of  Far  Outliers  approaches  10%,  the  analysis  should  be  repeated  on  logio  transformed  data. 
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It  should  not  be  assumed  a  priori  that  outliers  are  erroneous  data.  They  may  represent 
naturally  occurring  extreme  values,  e.g.  storm  runoff  concentrations. 

Beyond  identifying  extreme  values  for  closer  scrutiny,  outlier  elimination  on  successive  runs  is  intended  to 
improve: 

1.  graphics  scaling 

2.  non-robust  statistics,  particularly  the  estimates  of  percentage  variation  contributed  by  model  compo- 
nents. 

These  features  are  well  illustrated  in  the  example  discussed  below. 

5.3  EXAMPLES:  Logio  Transformation  and  Runs 

An  analysis  performed  with  program  TR2  of  a  particular  chloride  (CI)  data  set  is  illustrated  in  Figures 
11.1-11.4  with  no  data  transformation  and  in  Figures  11.5-11.8  for  logio  transformed  data.  Figures  11.1 
and  11.2  show  the  trend  plot  and  the  seasonal  plot  for  the  initial  Run  #1.  It  is  evident  in  both  cases  that 
the  overwhelming  majority  of  data  are  clustered  together  in  a  narrow  range  and  the  graphical  scale  has  been 
blown  up  substantially  to  meet  the  requirements  of  a  relatively  small  number  of  extreme  values.  From  the 
trend  plot  it  is  evident  that  virtually  all  of  these  have  occurred  since  1987  and  from  the  seasoned  plot  it  is 
clear  that  the  extremes  tend  to  occur  over  the  winter,  i.e.  late  in  the  year  and  early  in  the  following  year. 

Figures  11.3  and  11.4  show  the  comparable  plots  for  Run  #2  after  the  deletion  of  16  Fax  Outliers  that  were 
identified  during  the  initial  run.  The  scale  is  now  reduced  and  the  structure  in  both  trend  and  seasonal 
data  series  is  clearly  evident.  The  trend  plot  shows  an  abrupt  rise  in  CI  levels  that  appears  to  begin  in 
1987.  Although  rather  weak,  CI  seasonality  conforms  with  expected  behaviour  given  its  usage  as  winter  road 
de-icer.  A  review  of  the  16  outliers  in  the  report  file  shows  the  majority  of  them  to  have  occurred  during 
the  winter  months. 

The  other  feature  to  note  in  comparing  Run  #1  to  Run  #2  is  the  improvement  in  the  estimates  of  % 
Variability  contributed  by  trend  and  seasonal  components.  On  Run  #1  trend  contributed  only  about  9% 
of  the  variabihty  in  seasonally-adjusted  residuals  while  on  Run  #2  trend  contributed  about  51%  of  the 
variability.  Less  dramatically  on  the  seasonal  plots,  the  seasonal  function  contribution  to  variabihty  in 
de-trended  residuals  increases  from  about  3%  to  about  14%  after  deletion  of  the  16  Far  Outliers.  The  % 
Vciriability  estimates  are  analogues  of  the  R^  statistic  of  regression  einalysis  and  like  R^  aie  highly  sensitive 
to  the  presence  of  extreme  values.  One  of  the  approaches  to  developing  robust  R^  estimates  (e.g.  Rousseeuw 
aiid  Leroy,  1987)  has  been  to  weight  extreme  values  to  0  before  performing  the  calculations.  Zero  weighting 
is  virtually  equivalent  to  the  deletion  mechanism  that  is  built  into  the  TRENDS  software. 

Figure  11.5-11.8  show  the  comparable  analysis  with  program  TR2  performed  on  the  CI  data  series  after 
logio  transformation.  The  log  transform  seems  justified  here  as  the  data  clearly  span  nearly  2  orders  of 
magnitude.  Moreover,  the  results  obtained  after  the  initial  run  are  far  more  informative  than  those  obtained 
with  the  untransformed  series.  Run  #2  reveals  that  only  6  Far  Outliers  were  deleted,  another  indication 
that  the  logio  transform  is  the  better  choice.  The  mciin  improvement  of  Run  #2  occurs  with  the  seasonal 
function  where  scatter  and  consequently  the  ordinate  scale  has  been  reduced  by  fully  1  order  of  magnitude. 

5.4  ITERATIONS:  Seasonal-adjustment  models  (TRX  &:  TR2  only) 

The  need  for  an  iterative-decomposition  scheme  to  fit  seasonal  adjustment  models  arises  because  the  initial 
guesses  at  trend  term  T,  and  seasonal  term  5,  may  be  rather  crude.  Furthermore  the  problems  may  be 
aggravated  by  chronologically  or  seasonally  biased  data  distributions  at  various  times  in  the  record.  The 
benefits  of  iterative  fitting  are  illustrated  in  the  following  example  of  a  lake  total  phosphorus  (P)  series. 

5.4.1  ITERATIONS  EXAMPLE:  Heart  Lake  Total  P,  1975-1988. 

This  example  concerns  a  nutrient  data  series  gathered  in  a  long  term  lake  eutrophication  study.  From 
1975-1988  samples  were  collected  almost  exclusively  from  May-October  during  the  study  years.  Almost  no 
samples  were  collected  during  the  November-April  period.  Due  to  the  regularly  recurrent  5  month  gap,  the 
program  TRX  is  considered  the  most  appropriate  for  analyzing  the  historical  P  records.   An  initial  round 
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of  analysis  identified  2  extraordinarily  high  data  points  that  were  deleted  to  produce  a  simulation  of  what 
TRX  would  produce  during  the  actual  Run  #2  over  successive  iterations. 

The  effects  of  iterative  fitting  the  trend  component  are  illustrated  on  the  sequence  of  Figures  12.1-12.5 
which  correspond  respectively  to  the  output  after  1,  2,  3,  5  and  10  iterations.  Since,  on  the  first  iteration 
the  seasonal  term  is  effectively  assumed  to  be  zero,  Figure  12.1  shows  trend  determined  on  the  raw  data.  As 
we  see,  from  Figure  12.1  to  12.2,  there  is  a  drajnatic  reconfiguration  of  the  trend  line  on  the  second  iteration 
where  the  data  have  been  adjusted  by  the  first  guess  at  the  seasonal  function.  Also  note  that  the  estimates 
of: 

•  the  mean  of  the  annual  medians  has  dropped  somewhat  from  55.2  to  51.6  ^lg/L 

•  the  maximum  of  the  annual  medians  has  dropped  significantly  from  90  to  77.7  /ig/L 

•  the  estimate  of  %  variability  contributed  by  the  trend  term  to  the  de-seasonalized  residuals  has  virtually 
doubled  from  10.8  to  21.1%. 

In  going  from  Figure  12.2  to  12.3,  and  finally  from  Figure  12.3  to  12.4  further  adjustment  is  observed.  A 
final  run  with  20  iterations  produced  no  noticeable  change  over  the  results  obtained  at  10  iterations  which 
differed  only  slightly  from  results  obtained  with  5  iterations.  A  review  of  the  report  files  shows  that  in  4 
years  out  of  14,  the  final  secisonally  adjusted  annual  median  figures  after  10  iterations  declined  by  more  than 
10  f^g/L  from  the  original  estimates  of  annual  median  based  on  un-adjusted  data. 

Figures  12.6,  12.7  £ind  12.8  show  the  seasonal  plots  corresponding  to  the  trend  plots  of  Figure  12.1,  12.2 
and  12.5.  It  is  clear  that  for  the  P  data,  the  seasonal  function  on  the  first  iteration  (Figure  12.6)  gives  a 
reasonably  good  approximation  of  the  seasonal  periodicity.  On  the  second  iteration,  the  msdn  difference  is 
that  the  seasonal  peak  has  shifted  forward  by  a  few  days  and  the  data  exhibit  less  scatter  about  the  smoothed 
seasonal  median  which  accounts  for  the  increase  in  the  %  variability  explained  by  the  seasonal  function.  From 
the  second  iteration  forward,  the  changes  in  seasonal  function  are  scarcely  perceptible;  although,  the  seasonal 
amplitude  increases  slightly  to  reach  its  final  figure  of  40.6  /ig/L  at  the  5th  iteration  at  which  the  the  % 
explained  variability  also  stabilizes  at  about  49.3%. 

The  changes  observed  over  successive  iterations,  particularly  in  the  trend  estimates,  owe  their  explanation 
to  sampling  biases  embedded  in  individual  year's  data  relative  to  the  net  secisonal  characteristics  exhibited 
by  the  data.  For  the  4  years  of  record  for  which  annual  median  estimates  declined  appreciably  on  iteration, 
the  av^lilable  observations  were  biased  towards  the  peak  the  annual  seasonal  cycle  and  were  consequently 
adjusted  downwards. 

Although  the  seasonally  adjusted  annual  median  medians  show  appreciable  shifting  as  iterations  progress,  it 
should  be  noted  that  there  remains  considerable  scatter  in  the  adjusted  data;  such  that,  rigorous  statistical 
testing  may  not  show  significant  year-to-year  differences  between  many  of  the  years  which  have  altered 
position  after  the  first  few  iterations. 


6.0  TRENDS  Secondary  Program  Control  Parameters: 

6.1  Data  Series  Bracketing: 

Line  5  of  all  programs  .CNT  files  allows  users  to  enter  a  pair  of  'bracketing  dates'  which  can  be  used  to 
delimit  a  segment  of  the  input  data  series  for  analysis.  Only  samples  falling  between  the  two  dates  are 
accepted  for  analysis.  Reasons  for  using  bracketing  dates  include: 

1.  excluding  older  records  of  dubious  character  (see  Example  1), 

2.  excluding  stray  points  far  removed  from  the  main  series  and  particularly  with  program  TRX,  excluding 
years  at  the  beginning  or  end  of  record  with  few  data  points  (see  Example  2), 

3.  restricting  length  of  very  large  data  series  to  meet  program  size  constrciints. 

Generally,  on  preliminary  runs  the  bracketing  dates  should  be  set  to  span  the  full  time  horizon  of  the  available 
data  records. 
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6.1.1  Bracketing  Example  1: 

Figures  13.1  and  13.2  illustrate  the  potential  use  of  bracketing  dates.  The  original  bracketing  dates  were 
730000  and  910000.  Figure  13.1  shows  a  Na  series  for  which  the  first  1  1/2  years  of  record  comprise 
measurements  of  dubious  character  due  to  limitations  of  the  early  laboratory  techniques  available  at  the 
time.  A  check  of  the  input  data  file  revealed  that  reasonably  reliable  Na  records  begin  during  July  1974. 
The  final  observation  of  the  series  occurs  in  early  1990  and  is  the  only  datum  for  that  year.  For  a  data  set 
comprising  only  a  single  point,  the  median  defaults  to  the  value  of  that  point.  On  Figure  13.1,  the  final 
Jinnual  median  for  1990  suggests  a  continuing  upward  trend  which  is  not  justifiable  on  the  basis  of  a  only 
a  single  observation.  Figure  13.2  was  produced  by  rerunning  the  analysis  with  the  bracketing  dates  740700 
and  900000. 

6.1.2  Bracketing  Example  2: 

Figures  14.1  and  14.2  illustrate  another  example  for  which  series  bracketing  is  used  with  program  TRX.  In 
Figure  14.1,  the  first  datum  occurs  during  1973,  but  the  record  begins  in  earnest  only  during  1975.  Also 
at  the  right  end  of  the  Fe  series,  the  final  datum  occurs  during  early  1991  and  is  the  only  datum  available 
for  that  year.  Hence  at  both  extremities  of  the  Fe  series,  program  TRX  determines  annual  medians  for 
1973  and  1991  based  only  on  1  datum  in  each  case.  As  Figure  14.1  reveals  this  can  create  the  misleading 
impression  that  Fe  levels  have  dropped  dramatically  in  1990.  Figure  14.2  produced  by  using  the  bracketing 
dates  7500000  and  900000  yields  a  more  reliable  visual  impression. 

6.2  Zero  Option  Codes: 

Line  6  of  all  programs  .CNT  files  permits  users  to  choose  from  4  ways  of  dealing  with  0  values  that  might 
be  present  in  the  data  set.  Except  for  a  few  water  quality  variables,  e.g.  temp>erature  and  Gran  alkalinity, 
O's  should  not  be  present  in  water  quality  data  sets.  Rather  than  0,  zero  measurements  should  be  reported 
as  <  some  quantitative  value,  e.g.  a  detection  limit.  Nevertheless  O's  may  sometimes  be  encountered  in  older 
records  and  data  sets  prepared  by  others.  If  there  axe  only  a  few  O's,  the  user's  choice  amongst  the  4  available 
options  may  have  little  effect  on  the  results.  Since  logarithms  cannot  be  taken  of  O's,  the  programs  will  stop 
and  issue  an  error  message  if  the  control  file  indicates  the  Logio  transform  and  O's  are  encountered. 

6.3  Missing  Value  Code: 

The  missing  value  code  [MVC]  of  Line  6  of  all  programs  permits  the  user  to  input  data  from  data  matrix 
files  where  missing  observations  have  been  coded  with  an  MVC  as  discussed  in  Input  Data  File  Formats. 

6.4  Multiplier/Divisor  Codes: 

The  Multiplier/Divisor  options  of  Line  7  permit  users  to  change  data  to  more  appropriate  units  for  analysis. 
The  following  two  common  cases  arise: 

1.  awkward  data  retrieved  at  SIS  default  specifications 

2.  conversion  of  ion  data  from  mass  concentrations  to  charge  equivalents,  i.e.  fig/L  to  /^eq/L 

6.4.1  Multiplier/Divisor  Example  1:  Simple  unit  conversions 

Case  1  is  occurs  most  commonly  when  users  neglect  to  specify  units  when  retrieving  data  from  SIS.  Perhaps 
the  most  frequent  case  arises  with  the  variable  total  phosphorus  (P)  which  SIS  provides  at  default  units  of 
mg/L.  Except  for  heavily  polluted  sites,  P  levels  in  rivers  and  lakes  are  generally  in  the  /ig/L  range.  By 
working  in  mg/L  units  the  analyst  must  carry  additional  O's,  a  decimal  point,  and  if  taking  logarithms, 
negative  values  all  of  which  contribute  visual  clutter  on  the  y-axis  of  graphics  as  illustrated  in  Figure  15.1 
(on  the  Figure  15.1  y-cixis,  the  particular  scaling  factor  '♦lO"''  is  an  uncontrollable  foible  of  the  PLOT88 
plotting  software).  By  employing  the  Multiplier/Divisor  options  to  convert  the  input  P  series  to  /ig/L  units, 
Figure  15.2  has  much  cleaner  y-ajcis  labelling  and  annotation. 
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Figures  15.3  and  15.4  illustrates  similar  effects  of  unit  selection  on  certain  logarithmic  plots.  When  Fe  is 
expressed  in  mg/L,  the  j/-axis  tic  mark  labels  are  forced  to  carry  the'-'  sign  as  opposed  to  the  more  natural 
expression  in  /ig/L  on  Figure  15.4. 

6.4.2  Multiplier/Divisor  Example  2:  Conversion  to  charge  equivalents 

Case  2  probably  would  arise  most  frequently  when  dealing  with  the  dilute  waters  of  the  Canadian  Shield  where 
knowledge  of  fundamental  ionic  strength  indicators  and  the  principle  of  electroneutreility  are  frequently  of 
critical  importance.  Since  most  MOE  ion  data  are  reported  in  mass  concentration  units  and  the  SIS  retrieval 
options  do  not  allow  for  conversion  to  charge  equivalents,  the  Multiplier/Divisor  options  can  be  employed  to 
convert  units.  Figure  16.1  was  generated  using  the  sulphate  (SO^")  data  as  supplied  in  mass  concentration 
units.  Developing  the  results  in  charge  equivalent  units  may  be  done  as  follows.  First  obtain  the  equivalent 
weight  of  S04~  =  96.0576/2  =  48.0288.  Supposing  that  we  wish  to  obtain  the  results  in  peq/L  from  mg/L, 
we  get  the  number  of /ieq/L  per  mg/L  as  1000/48.0288  =  20.82084.  Then  the  Multiplier  option  is  invoked  by 
placing  'M'  in  column  1  of  Line  7  of  the  .CNT  file  and  multipher  20.82084  in  columns  2-10.  Equivalently  we 
could  proceed  by  invoking  the  Divisor  option  by  placing  'D'  in  column  1  and  the  divisor  .0480288  in  columns 
2-10.  The  results  are  shown  in  Figure  16.2  where  the  '*10^'  scaling  factor  on  the  y-axis  is  an  uncontrollable 
feature  of  the  PLOT88  plotting  software  scaling  routines.  Of  course,  if  a  data  series  were  supplied  in  charge 
equivalent  units  and  mass  concentration  units  were  desired,  the  Multiplier/Divisor  options  could  be  invoked 
perform  the  conversion. 

6.5  Plot  Codes: 

Plot  codes  are  entered  in  the  .CNT  files  on  Line  10  programs  (TRX  k  TRl)  and  Line  11  program  TR2. 
As  indicated  earlier,  usage  of  the  first  plot  code  'P'  to  generate  laser  printer  graphics  is  straight  forward. 

6.5.1  Plot  Density  Code 

The  plot  density  code  entered  in  column  4  allows  users  to  control  the  quality  of  the  graphics  produced.  Usage 
is  illustrated  in  Figures  17.1-17.3.  Figure  17.1  is  a  single  density  plot  obtained  by  placing  'S'  in  column  4. 
Figure  17.2  is  a  double  density  plot  obtained  by  placing  'D'  in  column  4.  Figure  17.3  is  a  quad  density 
plot  obtained  by  placing  'Q'  in  column  4.  The  differences  in  plot  quality  are  readily  evident.  It  is  important 
to  note  that  the  increased  quality  of  higher  density  plots  is  rather  inversely  proportional  to  the  time  required 
to  generate  the  plot.  Regardless  of  the  power  of  the  computer  used  the  speed  of  data  transmission  to  the 
laser  printer  remains  constant.  A  quad  density  plot  requires  the  transmission  of  approximately  4  times  the 
information  of  a  single  density  plot;  hence,  considerable  waiting  time.  To  reiterate  earlier  recommendations: 

•  Single  density  plots  are  appropriate  for  quick  rough  work.  Numerals  and  letters  may  be  rather  difficult 
to  decipher. 

•  Double  density  plots  are  acceptable  for  most  reporting  purposes.  Some  characters,  e.g.  '8',  are  not  well 
formed  and  may  suffer  further  on  reproduction.  Almost  all  graphics  in  this  manual  were  produced  at 
double  density. 

•  Quad  density  plots  are  recommended  for  high  quality  reporting  and  journal  article  graphics. 

6.5.2  Plot  Type  Code 

The  plot  type  code  entered  in  column  6  is  useful  for  controlling  the  type  of  trend  plot  generated  when  the  log 
transform  is  invoked.  When  no  log  transform  is  requested  the  plot  type  code  should  be  set  to  '0'  (an  error 
message  is  invoked  and  the  program  stopped  if  this  occurs).  Figures  16-17  are  examples  of  the  standard 
trend  plot  format  for  untransformed  data.  When  invoking  the  log  transform,  users  have  the  choice  of  the 
following  2  trend  plot  formats: 

1.  Set  to  '0'  a  normal  plot  is  produced  with  the  log  transformed  data  as  on  Figure  18.1.  Note  that  the 
j/-axis  is  scaled  identically  as  a  normal  plot. 

2.  Set  to  '1'  the  j/-axis  is  scaled  in  decadal  logio  cycles  as  on  Figure  18.2. 

The  reasons  for  offering  the  two  log  plot  types  should  be  evident  from  a  comparison  of  Figures  18.1  and 
18.2.   When  a  given  data  set  extends  across  two  or  more  decades  as  in  Figure  18.2,  it  may  occur  that  the 
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data  series  is  concentrated  into  the  centre  1/3  of  the  plot  leaving  about  2/3  of  the  plot  space  blank.  In 
contrast,  on  normally  scaled  Figure  18.1  the  data  series  fills  the  available  space  zoid  provides  a  superior 
visual  representation  of  the  series  behaviour  over  time.  Unfortunately,  due  to  the  erratic  scaling  routines  of 
the  PLOT88  software,  the  normally  scaled  log  plots  do  not  always  yield  such  satisfactory  results. 

Though  some  may  argue  for  exclusive  usage  of  decadal  style  plots  as  they  permit  direct  reading  of  numerical 
values,  illustrating  the  relative  shifts  in  the  data  series  over  time  is  the  primary  purpose  of 
these  graphics.  The  report  file  provides  summaries  of  annual  seasonally  adjusted  statistics. 

Microbiological  indicators,  e.g.  fecal  coliforms,  are  one  class  of  variables  which  typically  exhibit  such  high 
veiriability  that  the  decadal  style  log  plot  is  almost  always  most  appropriate,  as  in  Figure  18.3  on  which  the 
data  are  seen  to  extend  across  6  decades. 

6.5.3  Plot  Annotation  Code 

The  plot  annotation  code  in  column  8  allows  the  user  to  retain  or  delete  the  textual  information  and 
summary  statistics  appearing  on  the  right  side  of  the  plot.  An  'L'  placed  in  column  indicates  a  'long'  plot 
which  presents  the  information  on  the  right  side.  The  primary  reason  for  the  right  side  annotation  is  to  assist 
users  in  keeping  track  of  the  myriad  of  graphics  that  can  be  generated  rather  quickly  with  the  softweire. 

An  'S'  placed  in  column  indicates  a  "short"  plot  which  presents  only  the  plot  as  in  Figure  19.  The  primary 
reason  for  offering  this  option  is  to  facilitate  the  cutting  and  pasting  of  graphics  that  is  sometimes  required 
when  producing  reports  and  articles.  Note  that  the  two  plot  headers  and  the  y-axis  title  may  be  left  blank 
in  the  .CNT  file  as  discussed  below. 

6.5.4  Plot  Headers,  Run  Date  &:  1  -axis  Title: 

The  locations  of  the  2  plot  headers  (PLOT  HEADER  1  k.  PLOT  HEADER  2),  the  current  date  (TODAYS 
DATE)  and  the  y-axis  title  (Y-AXIS  TITLE)  are  illustrated  on  Figure  20.1  which  presents  the  TR2  analysis 
of  the  same  microbiologic£il  Figure  18.3.  Although  space  is  provided  for  these  features,  all  may  be  left  blank 
by  inserting  blank  lines  at  the  appropriate  locations  in  the  .CNT  files.  An  example  is  shown  in  Figure  20.3. 

6.6  Output  Files: 

Entering  any  non-'O"  digit  in  column  1  of  Line  16  (TRX  &:  TRl);  Line  17  (TR2)  of  the  .CNT  files 
generates: 

•  a  'Trend  Output  File'  for  all  3  programs 

•  a  'Seasonal  Output  File'  for  TRX  t  TR2 

according  to  the  names  k.  locations  specified  in  lines  17-18  of  the  .CNT  files.  The  features  of  these  output 
files  are  described  below. 

6.6.1  Trend  Output  Files: 
Line  1: 

•  The  first  18  columns  of  Line  1  of  trend  output  files  are  identical  to  those  of  the  input  data  file,  i.e.  the 
parameter  identifier  in  columns  1-6  and  the  site  identifier  in  columns  7-18. 

•  To  the  right  of  the  parameter  k.  site  identifiers,  beginning  in  column  23  appear  banners  identifying  the 
program  generating  the  file,  i.e. 

TRX  TREND  OUTPUT  FILE 
TRl  TREND  OUTPUT  FILE 
TR2  TREND  OUTPUT  FILE 

Lines  2— end  of  file:  Trend  file  output 

N.B.  All  trend  output  files  contain  the  same  number  of  lines  the  input  data  file  less  however  many  'far 
outliers'  have  been  deleted  by  the  final  'run'. 
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The  formats  of  trend  file  output  are: 

•  (I7,I5,F12.3,3F12.4)  program  TRl 

•  (I7,I5,F12.3,5F12.4)  programs  TRX  t  TR2 

The  1st  three  numerical  columns  of  all  trend  output  files  contain: 

•  sample  date  (yr/mo/dy)  in  columns  2-7 

•  sample  time  (hr/min)  in  columns  9-12 

•  simiple  concentration  A',  in  columns  13-24  in  format  F12.3 

The  4th  column  Jilways  contains  the  corresf)onding  trend  term  Ti  or,  if  the  analysis  has  been  performed  on 

logio  transformed  data,  antilog  T,  =  10  '. 

The  content  of  the  remaining  columns  to  the  right  vary  according  to: 

•  the  generating  program,  i.e.  TRl  output  differs  from  that  of  TRX  k.  TR2 

•  whether  the  analysis  was  performed  on  natural  or  log  transformed  input  data 

TRl  OUTPUT  TREND  FILE:  untransformed  data 

•  output  file  contains  only  4  numerical  columns  corresponding  to  the  sample  date,  time,  concentration  A', 
and  the  trend  estimate  T, 

TRl  OUTPUT  TREND  FILE:  logio  transformed  data 

•  the  1st  4  numerical  columns  contain  the  sample  date,  time,  concentration  A",  and  the  antilogged  trend 
estimate  10-^' 

•  numerical  column  5  contains  log,o  A', 

•  numerical  column  6  contains  T,  the  trend  estimate  determined  on  logjg  A\ 

TRX  &:  TR2  OUTPUT  TREND  FILE:  untransformed  data 

•  the  1st  4  numerical  columns  contains  the  sample  date,  time,  concentration  A',  and  the  trend  estimate 
T, 

•  numerical  column  5  contains  S,  the  relative  seasonal 

•  numerical  column  6  contains  T,  +  5,,  i.e.  the  sum  of  the  systematic  components  estimated  by  the  time 
series  model 

TRX  &  TR2  OUTPUT  TREND  FILE:  logic  transformed  data 

•  the  1st  4  numerical  columns  contain  the  sample  date,  time,  concentration  A',  and  the  antilogged  trend 
estimate  10^' 

•  numerical  column  5  contains  lO''^'''''^' ,  i.e.  the  antilogged  sum  of  the  systematic  components  estimated 
by  the  time  series  model 

•  numerical  column  6  contains  logjg  A", 

•  numerical  column  7  contains  T,  the  trend  estimate  determined  on  logjg  A', 

•  numerical  column  8  contains  S,  the  seasonal  estimate  determined  on  logjo  A'j 

6.6.2  Trend  Output  File  Use: 

The  trend  output  files  have  several  possible  uses  including: 

•  detailed  examination  of  the  numerical  values  of  T,,  Si  and  other  comjwnents  for  specific  samples  or  time 
periods 

•  selective  graphics  regeneration 

•  sources  files  for  compilation  of  special  graphics  in  other  software  packages,  e.g. 

correcting  PLOT88  idiosyncratic  scaling 
same  site,  multi-variable  overlays 
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sajne  variable,  multi-site  overlays 

6.6.2.1  Selective  Graphics  Regeneration: 

Presume  that  the  preliminary  analysis  of  a  water  quality  record  is  allowed  to  run  to  completion,  i.e.  until  no 
further  far  outliers  are  deleted,  and  quick  single  density  graphics  are  generated.  Suppose  that  completion  is 
achieved  on  the  4th  'run'.  To  later  regenerate  report  quality  quad  density  graphics  for  the  4th  run  without 
rerunning  all  graphics  the  following  procedure  may  be  applied: 

1.  if  there  exists  a  trend  output  file  from  the  first  analysis: 

-  specify  this  as  the  input  data  file 

-  run  the  desired  program  requesting  quad  density  plots 

2.  if  there  does  not  exist  a  trend  output  file  from  the  first  analysis: 

-  rerun  analysis  specifying  4  runs,  no  graphics  and  output  files 

-  specify  the  trend  output  file  as  the  input  data  file 

-  run  the  desired  program  requesting  quad  density  plots 

These  procedures  can  save  considerable  time  and  paper  when  selected  high  queility  graphics  must  be  regen- 
erated for  reports  or  journal  articles. 

6.6.3  Seasonal  Output  Files: 

The  'seasonal  output  files'  are  generated  only  by  programs  TRX  and  TR2. 
Line  1: 

•  The  first  18  columns  of  Line  1  of  seasonal  output  files  are  identical  to  those  of  the  input  data  file,  i.e.  the 
parameter  identifier  in  columns  1-6  and  the  site  identifier  in  columns  7-18. 

•  To  the  right  of  the  parameter  k  site  identifiers,  beginning  in  column  23  appear  banners  identifying  the 
program  generating  the  file,  i.e. 

TRX  SEASONAL  OUTPUT  FILE 
TR2  SEASONAL  OUTPUT  FILE 

Lines  2— end  of  file:  Seasonal  output  file 

The  secisonal  output  file  comprise  101  lines  with  2  numerical  columns  containing: 

1.  the  seasonal  decimal  date  series  from  0.00  to  1.00  year  in  steps  of  .01  year  at  which  the  seasonal  smooth 
was  determined 

2.  the  smoothed  relative  seasonal  function  values  for  the  respective  seasonal  decimal  dates 

The  seasonal  output  file  format  is;        (F7.2,F14.6) 

The  seasonal  output  files  are  useful  for: 

•  detailed  examination  of  the  seasonal  function 

•  source  files  for: 

-  comparative  analysis  of  sites  or  variables  via  graphic  overlays,  generally  within  other  software 
packages 

-  grouping,  e.g.  cluster  analysis,  of  multi-variable  or  multi-site  seasonal  functions 
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PART  III:  REPORT  FILES: 

1.0  Introduction: 

All  3  TRENDS  programs  generate  report  files  that  are  similar  in  structure  and  content,  but  vary  somewhat 
by  program  and  whether  the  logio  transform  has  been  invoked.  A  complete  set  of  6  report  files: 

•  TRX:  untransformed  series  • 

•  TRX:  log  10  transformed  series 

•  TRl:  untransformed  series 

•  TRl:  logio  transformed  series 

•  TR2:  untransformed  series 

•  TR2:  logio  transformed  series 

with  the  analysis  of  the  same  nitrate  series  is  available  in  the  Appendix. 

N.B.  In  the  .CNT  files  users  are  given  the  option  of  specifying  'short'  or  'long'  reports.  Generally  'short' 
reports  are  recommended.  The  'long'  report  includes  a  listing  of  input  data  and  correlations  amongst  model 
components.  Computation  of  the  latter  can  be  extremely  time  consuming.  The  correlations  are  currently 
under  investigation  as  a  tool  for  diagnosing  certain  apparently  anomalous  effects  that  can  be  observed  in  the 
approximate  %  contributions  to  total  series  variabihty  by  model  components  as  discussed  in  Section  6.2.4. 

2.0  General  Report  File  Layout: 

The  content  of  report  files  varies  according  to  the  time  series  model,  i.e.  by  program,  and  whether  or  not 
the  logio  transformation  is  invoked.  The  general  layout  is  as  follows: 

1.  Contents  of  the  program  control  (.CNT)  file  are  written  plus  additional  information  concerning  the 
number  of  observations.  Zero  observations  are  identified  and  their  disposition  noted.  If  invoked,  logio 
transforms  are  indicated.  The  Run  number  is  given.  For  programs  TRl  and  TR2,  a  summary  of 
sampling  interval  statistics  are  given  and  the  calculated  minimum  nearest  neighbours  to  be  employed 
with  the  smoothing  filter. 

2.  Trend  summary 

3.  Seasonal  summary 

4.  Approximate  %  contribution  to  total  variation  of  model  components  plus  summjiry  statistics 

5.  Model  residual  summary  statistics  and  outlier  identification  scan. 

When  the  'run'  parameter  is  >  1  and  far  outliers  have  been  identified,  deleted  and  the  trend  analysis  repeated, 
the  layout  of  the  following  'run'  is  identical  to  items  2-5. 

3.0  Definitions: 

The  following  definitions  are  critical  to  understanding  certain  of  the  summary  statistics  appearing  in  report 
files. 

N.B.  Note  the  following  significant  notational  change  from  the  conventions  employed  in  previous  sections: 

•  the  relative  seasonal  function  is  now  indicated  as  S,'  in  chronologic  time  or  SJ  in  seasonal  time  as  defined 
in  Section  7.2.3  below. 
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3.1  Time- Weighted  Mean  Trend: 

For  programs  TRl  and  TR2  employing  quasi-continuous  trend  models,  trend  summary  statistics  are  de- 
velof>ed  using  time— weighted  mean  trends  as  a  measure  of  central  tendency.  For  example,  the  grand 
time— weighted  mean  trend  T  for  the  series  is  determined  as 

2(tn-tl) 

where  T,  is  the  trend  estimate  for  chronologic  time  <,.  Similar  calculations  are  performed  to  obtiiin  the  annual 
time— weighted  mean  trends  for  the  year  by  year  trend  summaries.  The  antilogged  time-weighted  means 
trends  determined  on  logjo  transformed  data  are  referred  to  as  geometric  time— weighted  mean  trends. 

3.2  Relative  Trend: 

Use  of  the  relative  trend  series,  T,'  =  T,  —  T,  where  T  is  a  grand  measure  of  central  tendency  in  the  trend 
estimates  T,  over  the  time  horizon  of  the  series,  is  convenient  for  certain  algorithmic  manipulations  and  for 
the  determination  of  the  approximate  %  vciriation  contributed  by  the  trend  comp>onent.  As  indicated  for  the 
3  programs.  The  central  tendency  measure  is  taken  to  be: 

•  TRX:  T  =  mean  of  annual  medians  determined  on  seasonally  adjusted  residuals  r;,  =  A',  —  S/  where 
S,'  now  indicates  the  relative  seasonal 

•  TRl:  T  =  grand  time-weighted  mean  trend  T  determined  as  per  3.1  on  raw  or  transformed  series  A', 

•  TR2:  T  =  grand  time-weighted  mean  trend  T  determined  as  per  3.1  on  seasonally  adjusted  residuals 
t;,  =  Xi  —  S,'  where  S,'  now  indicates  the  relative  seasonal. 

3.3  Mean  Annual  Seasonal  Median: 

In  programs  TRX  and  TR2,  the  relative  seasonal  component  5/  is  not  determined  exactly  as  was  indicated 
earlier  in  PART  I:  Section  4.2.1,  but  rather  by  the  following  computationally  equivalent  approach: 

1.  obtain  de-trended  residuals  u,  —  X,  —  T/  where  7^'  is  the  relative  trend 

2.  reorder  de-trended  residuals  i/,  by  time  within  the  year  to  obtain  the  series  i/j 

3.  apply  MSM  filter  to  Vj  to  obtain  trial  seasonal  Sj 

4.  obtain  mean  annual  secisonal  median  as  5  =  5Zi=i  «5j/100 

5.  obtain  relative  seasonal  as  S,'  =  Sj  —  S 

6.  reorder  relative  seasonal  S'  to  chronologic  relative  seasonal  S,' 

This  approach  yields  the  mean  annual  seasonal  median  estimate  S  which  is  a  grand  central  tendency  measure 
for  the  series  that  is  free  of  seasonal  sampling  bias.  Generally  after  a  few  iterations  S  is  close  to  the  measure 
T  described  above.  Though  convergence  is  extremely  slow,  if  allowed  to  continue  for  many  iterations,  S  and 
T  graidually  become  identical  for  all  practical  purposes. 

4.0  Trend  Summaries: 

4.1  Program  TRX:  untransformed  data 

As  shown  on  the  example  in  the  Appendix,  for  untransformed  data  the  basic  TRX  trend  summary  comprises 
the  8  numerical  columns: 

1.  year  counter 

2.  year 

3.  number  of  observations  for  that  year 


PART  III:  REP  OR  T  FILES  32 

4.  annual  median  of  seasonally  adjusted  data,  iteration  =  1 

-  N.B.  recall  from  PART  I:  4.2.1  k.  PART  II:  5.4,  that  on  iteration  1  the  seasonal  is  assumed  to  be 
zero;  hence,  column  4  contains  the  annual  medians  of  raw  data 

5.  annual  median  of  seasonally  adjusted  data,  iteration  =  ITER,  the  number  of  iterations  specified  by  the 
user  in  the  .CNT  file 

6.  the  difference  between  the  annual  median  at  iteration  =  1  and  the  annual  median  at  iteration  =  ITER 

-  users  can  examine  the  list  to  determine  the  effects  of  iterative  decomposition 

-  large  differences  in  any  year  are  indicative  of  significant  seasonal  sampling  bias  for  that  yea^ 

7.  this  column  labelled  'YR-TO-YR\  simply  indicates  the  year  to  year  difference  in  final  annual  medians 
that  appears  in  column  8 

8.  the  difference  between  annual  medians  for  the  years  indicated  in  column  7 

-  provided  for  users  convenience  in  determining  successive  year  to  year  jumps  of  the  trend  line 

Beneath  the  trend  summary  are  listed: 

•  the  means  of  the  annual  medians  for  iterations  1  and  ITER 

•  for  iteration  =  ITER,  the  maximum  annual  median,  the  minimum  annual  median  and  their  difference 

4.1.2  Spearman's  Rho  Trend  Test 

The  final  component  of  the  trend  summary  is  a  simple  Spearman's  Rho  trend  test  on  the  final  sequence 
of  annual  medians  obtained  on  iteration  =  ITER.  The  indicated  p(RHO)  value  which  also  appears  on  the 
graphics  is  the  significance  level  for  a  two  sided  test.  As  per  the  usual  conventions,  for  p/e.05,  Rho  is 
considered  to  be  significant  at  95%  confidence.  This  is  an  extremely  crude  test  that  should  not  be  relied  on 
alone.  For  programs  TRl  and  TR2,  the  test  results  should  only  be  considered  when  USI  =  1  and  possibly 
USI  =  2.  Future  editions  of  the  software  will  provide  additional  formal  tests. 

4.2  Program  TRX:  Logio  transformed  data 

The  TRX  summary  for  logio  transformed  data  provides  a  year  by  year  summary  of  the  results  as  determined 
in  the  logio  spcice  followed  by  a  similar  summary  of  the  antilogged  results.  The  cintilogged  mean  of  the 
annual  medians  determined  on  the  logarithms  is  the  geometric  mean  annual  median. 

4.3  Program  TRl:  untransformed  data 

The  TRl  trend  summary  comprises: 

1.  grand  time- weighted  mean  determined  over  the  quasi-continuous  trend  line 

2.  a  listing  of  successive  local  maixima  and  minima  over  the  trend  line  followed  by  the  overall  extrema  and 
their  date  of  occurrence 

3.  a  year  by  year  summary  of  annual  time- weighted  mean  trend  levels  similar  to  the  summary  of  program 
TRX 

4.  Spearman's  Rho  trend  test  of  the  sequence  of  annual  time- weighted  mean  trends 

4.4  Program  TRl:   Logio  transformed  data 

Similar  to  untransformed  summary,  but  with  important  results  replicated  after  antilogging. 

4.5  Program  TR2:  untransformed  data 

Like  the  TRl  summary  the  TR2  summary  presents: 

1.  grand  time-weighted  mean  determined  over  the  quasi-continuous  trend  line 

2.  a  listing  of  successive  local  maxima  and  minima  over  the  trend  line  followed  by  the  overall  extrema  and 
their  date  of  occurrence 

3.  a  year  by  year  summary  of  annual  time-weighted  mean  trend  levels 

4.  Speajman's  Rho  trend  test  of  the  sequence  of  annual  time-weighted  meam  trends 


PART  III:  REPORT  FILES  33 

4.6  Program  TR2:   Logio  transformed  data 

Similar  to  untransformed  summary,  but  with  important  results  replicated  after  antilogging. 

5.0  Seasonal  Summaries: 

Seasonal  summaries  appear  only  in  report  files  generated  by  programs  TRX  and  TR2.  Since  the  same 
procedures  are  employed,  the  formats  are  identical. 

5.1  Seasonal  Summaries:  untransformed  data 

The  summary  comprises; 

1.  the  mean  annual  seasonal  median  as  per  Section  3.3 

2.  a  scan  of  local  extrema  over  the  seasonal  trend  line 

3.  summary  of  global  extrema  over  the  seasonal  trend  line  including  the  seasonal  amplitude  determined  as 
the  annual  maximum  —  annual  minimum 

5.2  Seasonal  Summaries:  Logic  transformed  data 

The  summary  comprises; 

1.  the  mean  annual  seasonal  median  as  per  Section  3.3  in  logio  units  and  backtransformed  to  natural  units 

2.  a  scan  of  local  extrema  over  the  seasonal  trend  line 

3.  summaries  of  global  extrema  over  the  sefisonal  trend  line  including  the  seasonal  amplitude  determined 
as  the  annual  maximum  —  annual  minimum 

these  are  determined  at  three  levels: 

a.  geometric  mean  annual  seasonal  median 

b.  maximum  annual  time-weighted  mean  trend 

c.  minimum  annual  time-weighted  mean  trend 

so  that  analysts  may  examine  how  seasonal  amplitude  changes  with  trend  level  with  the  geometric 
model,  i.e.  logio  transformed  data 

6.0  Approximate  %  Contributions  to  Variation  &i  Summary  Statistics; 

The  summary  reports  for  al\  programs  contain  estimates  of  the  approximate  contributions  to  variation  by 
model  components.  These  estimates  reflect  variability  about  grand  robust  central  tendency  measures  of  the 
series  A'i,  hence  approximate.  Also  a  subtractive  rather  than  additive  principle  is  used. 

6.1  TRl 

For  program  TRl,  approximate  %  contributions  to  total  veo-iation  in  the  raw  or  trainsformed  data  series  A", 
are  determined  as  follows: 

•  total  variation  =  13"=  i(-'^'  ~  ^)^ 

•  variation  due  to  trend  =  ^"=i(A',  -  ff  -  E"=i(^'i  -  f  -  T,')^ 

•  vciriation  due  to  residual  noise  =  total  variation  —  variation  due  to  trend 

where  the  series  is  considered  to  be  centred  about  the  grand  time-weighted  mean  trend  T.  The  variation 
due  to  trend,  when  expressed  as  a  percentage  of  total  variation  is  equivalently  an  approximate  R^  measure 
for  the  model. 

6.2  TRX  &:  TR2 

For  programs  TRX  k.  TR2  approximate  contributions  to  variation  are  determined  for  the  3  series: 
1.  de-seasonalized  residuals  -qi  =  A',  —  5/ 
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-  relative  trend  component  T,' 

-  noise  component  e,  =  t],  —  T  —  T/ 

2.  de-trended  residuals  2/,  =  A';  —  Tj' 

-  relative  sectsonal  component  S,' 

-  noise  component  e,  =  r;,  —  S  —  S/ 

3.  raw  or  transformed  series  A\ 

-  relative  trend  component  T,' 

-  relative  seasonal  component  S/ 

-  relative  systematic  components  T,'  -I-  S/ 

-  noise  component  £;  =  A',  —  A'  —  T,'  —  S,' 

where: 

•  relative  trend  T/  and  central  tendency  measure  T  are  determined  as  in  Section  3.2 

•  relative  seasonal  5/  and  central  tendency  measure  S  are  determined  as  in  Section  3.3 

•  central  tendency  measure  A"  is  taken  as  {T  +  S)/2 

6.2.1  Variation  in  de-seasonalized  series: 

Approximate  %  contributions  to  total  variation  in  the  seasonally  adjusted  series  r],  =  A\  —  5/  are  determined 
as  follows: 

•  total  variation  =  Yl't  =  i(^>  ~  ^)^ 

•  variation  due  to  trend  =  X^"=,(r?,  -  T)^  -  IZ,"=i(''>  -T  -  T,')^ 

•  variation  due  to  residual  noise  =  total  variation  —  variation  due  to  trend 

6.2.2.  Variation  in  de-trended  series: 

Approximate  %  contributions  to  total  variation  in  the  seasonally  adjusted  series  i/,  =  Xi  —  T|'  are  determined 
as  follows: 

•  total  variation  =  Yl^=ii^>  ~  ^)^ 

•  variation  due  to  trend  =  Yl'i=ii^'  "  ^V  -  T,"=ii''i  "  ■^  -  S,')^ 

•  variation  due  to  residual  noise  =  total  variation  —  variation  due  to  secisonal 

6.2.3  Total  variation  in  raw  or  transformed  series: 

Let  X  =  {T  +  S)/2  be  a  grand  robust  central  tendency  measure  of  the  series  A',  over  its  time  horizon. 
Approximate  %  contributions  to  total  variation  in  the  raw  or  transformed  series  A',  are  determined  as 
follows: 

•  total  variation  =  ^"=i(A',  —  A')^ 

•  variation  due  to  trend  =  X^"^i(A',  -  X^  -  E"=i(-^''  -  A'  -  T,')^ 

•  variation  due  to  seasonal  =  13"^,(A',  -  A')^  -  X!"=i(A',  -  X  -  S,')^ 

•  variation  due  to  trend  +  seasonal  =  Y,'!=i{X,  -  X)^  -  ^,"^,(A',  -  X  -  T,'  -  5,')^ 

•  vairiation  due  to  residual  noise  =  total  variation  —  variation  due  to  trend  +  seasonal 

6.2.4  Interpreting  %  Contributions  to  Series  Variation: 

The  estimated  %  contributions  to  series  variability  are  approximate  figures  that  are  informative  but  should 
nevertheless  be  treated  with  caution.  Like  the  R'^  statistic  of  regression  analysis,  the  percentages  are  non- 
robust  and  highly  sensitive  to  the  presence  of  outliers.  Thus,  if  outliers  are  present,  the  %'s  only  begin  to 
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cissume  a  reasonable  degree  of  validity  after  2  or  3  "runs'  have  been  performed.  Also  recall  the  influence 
of  iterations  as  in  the  example  of  PART  II:  5.4.1  where  successive  iterations  significantly  increase  the  % 
contributions  to  variability. 

For  the  nitrate  series  of  the  sample  TRX  report  files  given  in  the  appendix,  readers  will  note  that  the  'TOTAL 
VARIATION  IN  SERIES"  summary  for  untransformed  data  reveals: 

•  variation  due  to  trend  =  58.3  % 

•  variation  due  to  seasonal  =  -5.8  % 

•  variation  due  to  trend  +  seasonal  =  68.4  %. 

Most  often  the  separate  trend  and  seasonal  contributions  should  add  to  produce  a  figure  approximately 
equal  to  the  joint  contribution.  Nevertheless,  the  above  figures  are  not  necessarily  in  error!  Rather  what 
they  reveal  are  certeun  synergisms  and  biases  amongst  components  embedded  in  the  series.  This  type  of 
situation  may  reflect: 

•  significant  correlation  between  estimated  trend  and  seasonal  components 

•  local  residual  variability  that  varies  over  time,  e.g.  the  case  of  older  records  that  show  much  more  scatter 
than  current  records 

•  local  residual  variability  that  varies  with  the  magnitude  of  the  underlying  process,  e.g.  variability  is 
positively  correlated  to  trend  level 

•  chronologically  varying  seasonal  sampling  biases,  e.g.  older  records  reflect  intensive  summer  sampling 
while  recent  records  are  evenly  distributed  throughout  the  year  when  a  strong  seasonality  is  present 

and  other  causes.  To  degrees  such  circumstances  represent  violations  of  the  general  assumptions  that  un- 
derlie the  time  series  models;  although,  since  the  models  are  inherently  robust,  results  concerning  trends 
do  not  suffer.  The  correlation  coefficients  that  are  determined  when  a  'long  report'  is  requested  are  un- 
der development  as  a  tool  for  diagnosing  the  causes  of  such  phenomena.  This  is  an  area  requiring  further 
research. 

6.3  Summary  Statistics: 

Just  beneath  the  approximate  %  contributions  to  variation  appear  summary  statistics  determined  on  the 
series  under  examination  and  the  corresponding  residual  noise.  The  mean  and  standard  deviation  require  no 
explanation,  but  the  robust  statistics  may  be  unfamiliar  to  some.  Consider  robustness  to  imply  resistance 
to  distribution  shape  and  the  presence  of  extreme  values.  Definitions  follow: 

median: 

•  2nd  quartile  q-^  or  50th  quantile  Qso;  a  robust  measure  of  central  tendency 
MAD: 

•  median  of  absolute  deviations  from  a  robust  central  tendency  measure  A' 

•  herein  let  A'  =  median  {A\},  hence  MAD  A'  =  medicin  {|A',  —  X\). 

•  MAD  is  a  robust  measure  of  dispersion  ss  2/3<t  the  standard  deviation  when  data  are  normally  dis- 
tributed 

IQR: 

•  interquartile  range,  i.e.  93  -  qi  where  93  =  3rd  quartile  =  Q75  =  75th  quantile  and  qi  =  1st  quartile  = 
Q25  —  25th  quantile 

•  IQR/2  is  a  robust  measure  of  dispersion  s;  2/Zcr  the  standard  deviation  when  data  are  normcdly  dis- 
tributed 

Robust  Standard  Deviation: 

•  a'  =  5'/$-'(.75)  -  1.48260225*  where  5*  is  a  robust  measure  of  dispersion  such  as  MAD  or  1/2  IQR 

•  the  reported  robust  standard  deviation  is  derived  with  MAD 
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•  N.B.  4'(-)  is  the  normal  distribution  function,  thus  a'  assumes  that  robust  dispersion  S'  is  derived  from 
data  that  over  the  central  50%  of  their  distribution  behave  approximately  normally 

N.B.  If  all  A',  are  of  the  same  sign,  median  {A',}  and  median  {|A',|}  will  be  identical 

7.0  RESIDUAL  SUMMARY  STATISTICS  ii  OUTLIER  IDENTIFICATION: 

Finally  for  the  model  residuals  £, ,  basic  statistics,  including  order  statistics,  i.e.  quantiles,  and  remge  statistics, 
are  summarized  and  an  outlier  analysis  is  performed  to  categorize  extreme  residuals  by  several  classifications 
in  common  use.  Far  outliers  identified  as  exceeding  Tukey's  (1977)  outer  fences  £ire  deleted  from  series  A', 
in  preparation  for  a  repkcat  trend  analysis  on  a  subsequent  'run'.  An  outer  fence  is  75'  where  S*  =  IQR/2 
which  for  normally  distributed  data  is  about  4.72  standard  deviations  from  the  mean. 
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APPENDIX  B 
TABLES 


for 
PC-Software,  Users  Guide  and  Documentation 

for 

Robust  Graphical  Time  Series  Analysis 

of 

Long  Term  Surface  Water  Quality  Records 


APPENDIX:  TABLE  1 


Table  1.  Robust  Smoothing  Filters  of  Programs  TRl  and  TR2: 

Filter  1: 

•  basic  LOWESS  with  2  robustness  iterations 

FUter  2: 

•  LOWESP  (modified  LOWESS)  with  2  robustness  iterations 

Filter  3: 

3a.  basic  LOWESS  with  2  robustness  iterations 

3b.  basic  LOWESS  with  0  robustness  iterations  on  output  of  3a. 

3c.  basic  LOWESS  with  0  robustness  iterations  on  output  of  3b. 

Filter  4: 

4a.  LOWESP  with  2  robustness  iterations 

4b.  LOWESP  with  0  robustness  iterations  on  output  of  4a. 

4c.  LOWESP  with  0  robustness  iterations  on  output  of  4b. 

Filter  5: 

5a.  LOWESP  with  2  robustness  iterations. 
5b.  output  of  5a  is  interpolated  to  an  evenly  spaced  series. 
5c.  LOWESP  with  2  robustness  iterations  on  evenly  spaced  series  5b. 
5d.  output  of  5c  is  interpolated  to  original  irregular  series. 
5c.  LOWESP  with  0  robustness  iterations  on  irregular  series  5d. 
Filter  6: 

6a.  LOWESP  with  2  robustness  iterations. 

6b.  LOWESP  with  0  robustness  iterations  on  output  6a  as  if  series  were  evenly  spaced. 
6c.  LOWESP  with  0  robustness  iterations  on  output  6b  at  original  irregular  sp>acing. 
Filter  7: 

7a.  Tukey  H3RSSH  smooth  as  if  series  were  evenly  spaced 
7b.  Tukey  H3RSSH  smooth  of  7a  output  as  if  series  were  evenly  spaced  . 
7c.  LOWESP  with  2  robustness  iterations  on  output  7b  at  original  irregiilar  spacing. 
7d.  LOWESP  with  0  robustness  iterations  on  output  7c  at  original  irregular  spacing. 
FUter  8: 

8a.  Tukey  H3RSSE  smooth  as  if  series  were  evenly  spaced 
8b.  Tukey  H3RSSH  smooth  of  8a  output  as  if  series  were  evenly  sp>aced 
Filter  9: 

•  basic  LOWESS  with  0  robustness  iterations 


APPENDIX:  TABLE  2 


Table  2.   Seasonal  Date  Conversion  Chart, 
(365  day  year) . 


Seasonal 

Caler 

idar 

Date 

Date 

Time 

.00 

JAN 

1 

00 

00 

.01 

JAN 

4 

15 

36 

.02 

JAN 

8 

07 

12 

.03 

JAN 

11 

22 

48 

.04 

JAN 

15 

14 

24 

.05 

JAN 

19 

06 

00 

.06 

JAN 

22 

21 

36 

.07 

JAN 

26 

13 

12 

.08 

JAN 

30 

04 

48 

.09 

FEB 

2 

20 

24 

.10 

FEB 

6 

12 

00 

.11 

FEB 

10 

03 

36 

.12 

FEB 

13 

19 

12 

.13 

FEB 

17 

10 

48 

.14 

FEB 

21 

02 

24 

.15 

FEB 

24 

18 

00 

.16 

FEB 

28 

09 

36 

.17 

MAR 

4 

01 

12 

.18 

MAR 

7 

16 

48 

.19 

MAR 

11 

08 

24 

.20 

MAR 

15 

00 

00 

.21 

MAR 

18 

15 

36 

.22 

MAR 

22 

07 

12 

.23 

MAR 

25 

22 

48 

.24 

MAR 

29 

14 

24 

.25 

APR 

2 

06 

00 

.26 

APR 

5 

21 

36 

.27 

APR 

9 

13 

12 

.28 

APR 

13 

04 

48 

.29 

APR 

16 

20 

24 

.30 

APR 

20 

12 

00 

.31 

APR 

24 

03 

36 

.32 

APR 

27 

19 

12 

.33 

MAY 

1 

10 

48 

.34 

MAY 

5 

02 

24 

.35 

MAY 

8 

18 

00 

.36 

MAY 

12 

09 

36 

.37 

MAY 

16 

01 

12 

.38 

MAY 

19 

16 

48 

.39 

MAY 

23 

08 

24 

APPENDIX:  TABLE  2 


Table  2.  Seasonal  Date  Conversion  Chart  Cont'd, 
(365  day  year) . 


Seasonal 

Calendar 

Date 

Date 

Time 

.40 

MAY  27 

00:00 

.41 

MAY  30 

15:36 

.42 

JUN   3 

07:12 

.43 

JUN   6 

22:48 

.44 

JUN  10 

14:24 

.45 

JUN  14 

06:00 

.46 

JUN  17 

21:36 

.47 

JUN  21 

13:12 

.48 

JUN  25 

04:48 

.49 

JUN  28 

20:24 

.50 

JUL   2 

12:00 

.51 

JUL   6 

03:36 

.52 

JUL   9 

19:12 

.53 

JUL  13 

10:48 

.54 

JUL  17 

02:24 

.55 

JUL  20 

18:00 

•  .56 

JUL  24 

09:36 

.57 

JUL  28 

01:12 

.58 

JUL  31 

16:48 

.59 

AUG   4 

08:24 

.60 

AUG   8 

00:00 

.61 

AUG  11 

15:36 

.62 

AUG  15 

07:12 

.63 

AUG  18 

22:48 

.64 

AUG  22 

14:24 

.65 

AUG  2  6 

06:00 

.66 

AUG  29 

21:36 

.67 

SEP   2 

13:12 

.68 

SEP   6 

04:48 

.69 

SEP   9 

20:24 

.70 

SEP  13 

12:00 

.71 

SEP  17 

03:36 

.72 

SEP  20 

19:12 

.73 

SEP  24 

10:48 

.74 

SEP  28 

02:24 

.75 

OCT   1 

18:00 

.76 

OCT   5 

09:36 

.77 

OCT   9 

01:12 

APPENDIX:  TABLE  2 


Table  2.  Seasonal  Date  Conversion  Chart  Cont'd, 
(365  day  year) . 


Seasonal 

Calendar 

Date 

Date 

Time 

.78 

OCT 

12 

16: 

48 

.79 

OCT 

16 

08; 

24 

.80 

OCT 

20 

00: 

00 

.81 

OCT 

23 

15 

36 

.82 

OCT 

27 

07 

12 

.83 

OCT 

30 

22 

48 

.84 

NOV 

3 

14 

24 

.85 

NOV 

7 

06 

00 

.86 

NOV 

10 

21 

36 

.87 

NOV 

14 

13 

12 

.88 

NOV 

18 

04 

48 

.89 

NOV 

21 

20 

24 

.90 

NOV 

25 

12 

00 

.91 

NOV 

29 

03 

36 

.92 

DEC 

2 

19 

12 

.93 

DEC 

6 

10 

48 

.94 

DEC 

10 

02 

24 

.95 

DEC 

13 

18 

.00 

.96 

DEC 

17 

09 

.36 

.97 

DEC 

21 

01 

:12 

.98 

DEC 

24 

16 

:48 

.99 

DEC 

28 

08 

:24 

1.00 

DEC 

31 

24 

:00 

APPENDIX:  TABLE  3 


Table  3.  Sample  input  data  file. 

+ 1 + 2 + 3 

As     170026013 

170026013  641117     0  320. 

170026013  650308     0  130. 

170026013  650323     0  200. 

170026013  650405     0  50. 

170026013  650412     0  200. 

170026013  650531     0  80. 

170026013  650617     0  300. 

170026013  650621     0  300. 

170026013  650708     0  1000. 

170026013  650720     0  1500. 

170026013  650803     0  1000. 

170026013  650817     0  5000. 

170026013  650818     0  5000. 

170026013  650908     0  2500. 

170026013  650916     0  4000. 

170026013  650921     0  6000. 

170026013  651021     0  250. 

170026013  660830  1555  2500. 

170026013  671002  1430  5.< 

170026013  671101  1540  224. 

170026013  680402  1125  40. 

170026013  680522  1830  90. 

170026013  680604     0  490. 

170026013  680703     0  60. 

170026013  680807     0  540. 

170026013  680907     0  880. 

170026013  681125  1600  290. 

170026013  690217  1545  80. 

170026013  690314  1530  60. 

170026013  690422  1525  80. 

170026013  690618  1415  142. 

170026013  690715  1430  180. 

170026013  690819  1440  1000. 

170026013  690923  1430  1520. 

170026013  691021  1410  470. 

170026013  691111  1550  560. 


APf>QBIX:  TABLE  4 


Table  4.  Sanple  input  data  fUe  (Mtrix). 
.1.--  — ---2--- -♦--   -3-   ------4— ------5— - 


IONS  CLEARWATER 
730628  1200 
730718  1200 
730802  1200 
730809  1200 
730814  1200 
730821  1200 
730828  1200 
730906  1200 
730912  1200 
730919  1200 
730925  1200 
731003  1200 
731010  1200 
731015  1200 
731023  1200 
731102  1200 
731106  1200 
740123  1200 
740219  1200 
740403  1200 
740515  1200 
740522  1200 
740528  1200 
740604  1200 
740612  1200 
740618  1200 
740625  1200 
740703  1200 
740710  1200 
740717  1200 
740723  1200 
740730  1200 
740806  1200 
740815  1200 
740820  1200 
740827  1200 
740904  1200 
740910  1200 
740918  1200 
740925  1200 
741002  1200 
741008  1200 
741015  1200 
741114  1200 


87.5000 
91.0000 
91.0000 
86.5000 
85.0000 
86.0000 
91.0000 
99.0000 
87.0000 
92.0000 
90.0000 
90.0000 
92.0000 
91.0000 
91.0000 
89.3195 
87.0000 
93.0000 
98.5000 
97.0000 
86.0000 
87.5000 
86.0000 
82.0000 
84.0000 
S4.0000 
84.0000 
84.0000 
85.5000 
85.5000 
84.0000 
81.6670 
85.0000 
84.6670 
87.0000 
86.0000 
87.5000 
87.5000 
85.0000 
88.5000 
88.0000 
87.5000 
86.5000 
85.5000 


62.5950 
39.4948 
49.7209 
78.8024 
49.7209 
49.7209 
49.7209 
39.4948 
47.9798 
35.1998 
49.7209 
44.3139 
49.7209 
39.4948 
49.7209 
49.7209 
78.8024 
46.8597 
46.5178 
47.4875 
41.8021 
41.6212 
41.5290 
41.4964 
41.6734 
42.1103 
43.1132 
44.8428 
46.4174 
47.5864 
48.1177 
48.3363 
48.3U7 
48.2714 
48.2371 
48.2061 
48.1900 
48.1830 
48.1626 
48.0769 
47.8517 
47.5405 
47.1130 
45.7974 


-999.0000 
-999.0000 
-999.0000 
-999.0000 
-999.0000 
-999.0000 
-999.0000 
-999.0000 
-999.0000 
-999.0000 
-999.0000 
-999.0000 
-999.0000 
-999.0000 
-999.0000 
-999.0000 
-999.0000 
-999.0000 
-999.0000 
-999.0000 
-999.0000 
-999.0000 
-999.0000 
-999.0000 
-999.0000 
-999.0000 
-999.0000 
-999.0000 
-999.0000 
-999.0000 
-999.0000 
-999.0000 
-999.0000 
-999.0000 
-999.0000 
-999.0000 
-999.0000 
-999.0000 
-999.0000 
-999.0000 
-999.0000 
-999.0000 
-999.0000 
-999.0000 


-999.0000 
-999.0000 
-999.0000 
-999.0000 
-999.0000 
-999.0000 
-999.0000 
-999.0000 
-999.0000 
-999.0000 
-999.0000 
-999.0000 
-999.0000 
-999.0000 
-999.0000 
-999.0000 
-999.0000 
-999.0000 
-999.0000 
-999.0000 
-999.0000 
-999.0000 
-999.0000 
-999.0000 
-999.0000 
-999.0000 
-999.0000 
-999.0000 
-999.0000 
-999.0000 
-999.0000 
-999.0000 
-999.0000 
-999.0000 
-999.0000 
-999.0000 
-999.0000 
-999.0000 
-999.0000 
-999.0000 
-999.0000 
-999.0000 
-999.0000 
-999.0000 


-999.0000 

-999.0000 

-999.0000 

-999.0000 

-999.0000 

-999.0000 

-999.0000 

-999.0000 

-999.0000 

-999.0000 

-999.0000 

-999.0000 

-999.0000 

-999.0000 

-999.0000 

-999.0000 

-999.0000 

-999.0000 

-999.0000 

-999.0000 

-999.0000 

-999.0000 

-999.0000 

-999.0000 

-999.0000 

-999.0000 

-999.0000 

67.4212 

63.0715 

60.8966 

67.4212 

65.2463 

58.7217 

65.2463 

63.0715 

58.7217 

52.1971 

55.4594 

73.9458 

73.9458 

71.7710 

67.4212 

65.2463 

58.7217 


IX:  TMU  5 


Table  5a.   Sample  .CNT  file  for  prograai  TRX 
♦ 1 ♦----2 ♦----3 


D:\TRX\CRUD 
D:\ffloira\as011.dat 
As 
(10X,I7,l5,f10.1) 

730000  900000 
1        -9.000 
1  1.000 

1      1     5 

13     8 

P  S  1  L 
Arsenic  ug/L  (<OL) 
Moira  Hwy  62 
APR  22  1991 
Concentration  ug/L 
SR 
1 

D:\TRX\CRUD1 
O:\TRX\CRUD2 


1:report  file  name 

2:input  data  file  name 

3:input  variable  identifier 

4: input  data  fonaat 

S:bracket  dates 

6:zero  option  t   aisaing  value  codes 

7:inultiplier/di visor  codes 

8:log  code,  #  runs,  #  iterations 

9:seasonal  saoothing  parameters 

10:plot  codes 

11:plot  header  1 

12:plot  header  2 

13:current  date 

14:y-axi8  title 

15: long/short  report  code 

16:output  f i le  code 

17:output  trend  f i le 

18:output  seasonal  file 


Table  5b.  Sanple  .CNT  file  for  program  TR1 
♦ 1 ♦ 2 ♦ 3 


D:\TR1\CRUD 
D:\clear\ion.dat 
S04 
(17,I5,132X,f12.4) 

730000  910000 
1   -9.000 
H  20.82084 
0  3 

5  1  1. 

P  S  1  L 
Sulphate  ueq/L 
Clearwater  L 
HAY  11  1991 
Concentration  ueq/L 
SR 
1 
D:\TR1\CRUD1 


1 :report  f i le  name 

2: input  data  file  name 

3:input  variable  identifier 

4:  input  data  fortnat 

Srbracket  dates 

6:zero  option  t   Missing  value  codes 

7:Multiplier/divisor  codes 

6: log  code,  #  rins 

9:filter,  sample  interval  type,  USI 
10:plot  codes 
11:plot  header  1 
12:plot  header  2 
13:current  date 
14:y-axis  title 
15: long/short  report  code 
16:output  file  code 
17:output  trend  file 


APPQSIX:  TABU  5 


T«ble  5c.   Sample  .CUT  file  for  proflrim  TR2 
....... ..1--..+... -2-. --♦... -3 


D:\TR2\CRII) 
D:\clear\cond.dat 
CONO 

(5X,I5,132X,f9.3) 
730000  910000 


-9.000 


0 
5  1 
13 
P  0 


3  5 

1. 
8 

0  L 
Conductance 
Clearwater  I 
jin  18  1991 
Concentration  uS/cm 
SR 
1 

O:\TR2\CRU01 
D:\TR2\CRUD2 


1:report  f i le  naaK 

2:input  data  file  naae 

3:input  variable  identifier 

4:  input  data  fonaat 

S:bracket  dates 

6:zero  option  t  ditting  value  code* 

7:injltiplier/divisor  codes 

8:  log  code,  «  r\r» 

9:f liter,  saaple  interval  type,  USI 
10: seasonal  saoothing  parameters 
11:plot  codes 
12:plot  header  1 
13:plot  header  2 
14:current  date 
15:y-axis  title 
16: long/short  report  code 
17:output  f i le  code 
18:output  trend  f i le 
19:output  seasonal  file 


APPENDIX  C 
REPORT  FILES 


for 

PC-Software,  Users  Guide  and  Documentation 

for 

Robust  Graphical  Time  Series  Analysis 

of 

Long  Term  Surface  Water  Quality  Records 


APPOBIX:  KPCtT  FIL£S:  TIX  untrmfoTMd  dit« 


PROGRAM  TRX 

REPORT  DATA  FILE  «  D:\0ILDO\TRX\CfiHO3 
INPUT  DATA  FILE  «  O:\0ILDO\DAT2\NO3.dat 
VARIABLE  >  N03 
(5X,i7,i5,f9.3) 
START  TEAR  «  720000   ENO  YEAR  «  900000 

ZERO  OPTION  CODE  >  1 

NUNERIC  HISSING  VALUE  CODE  >    -9.00000000 

MULTIPLIER/DIVISOR  CODE  >  1 

MULTIPLIER/DIVISOR  =      .00000000 

INPUT  DATA  ARE  LEFT  AS  IS 

LOGG  CODE  =   0 

RUNS  =   1 

ITERATIONS  =   5 

MINIMUN  N  /  UINOOU  *   13 

MINIMUM  SEASONAL  SMOOTHING  UINDOU  ■  1/   8  yr 

PLOT  COOES: 
PLOT  INDICATOR  >  0 
PLOT  DENSITY  1  «  D 
LOG  PLOT  TYPE  *   0 

PLOT  TYPE  CODE  «  L   (i.e.  LONG/SHORT) 
PLOT  HEADER  1  >  Nitrate  ug/L 
PLOT  HEADER  2  =  CLEARWATER  L 
RUN  DATE  HEADER  '   JUN  16  1991 
Y-AXIS  TITLE  '   Concentration  us/L 
REPORT  FILE  CODE  =  SR  (i.e.  SHORT  REPORT) 
OUTPUT  FILE  CODE  »  0 
NO  OUTPUT  FILES  REQUESTED 
FIRST  LINE  OF  INPUT  DATA  FILE: 

PARAMETER  =  N03      SITE  =  CLEARWATER  L 

QUALITY  DATA  FOR  PARM  «  N03     READ  t  STORED 
212  RECORDS  WERE  READ 

THE  FOLLOWING  3  INPUT  DATA  VERl   IDENITIFIEO  AS  ZERO'S: 
870901  1100      .000  ZERO 
880830  1030     .000  ZERO 
880926  1000      .000  ZERO 

SMALLEST  NON-ZERO  POSITIVE  VALUE  ■     4.00000 

ZERO  OPTION  CODE  '   1 

SETTING  O'S  TO  1/2  SMALLEST  POSITIVE  NON-ZERO  VALUE  «     2.00000 


BEGINNING  RUN  1 


TREND 

SUWIARY 

YEAR 

N 

ITER  =  1 

ITER  >  5 

OIFF 

YR-TO-YR 

DIFF 

73 

16 

101.500 

110.371 

-8.871 

74-73 

-40.552 

74 

21 

61.000 

69.819 

-8.819 

75-74 

-.764 

75 

10 

66.000 

69.055 

-3.055 

76-75 

-12.524 

76 

13 

46.000 

56.531 

-10.531 

77-76 

3.996 

77 

10 

63.500 

60.528 

2.972 

78-77 

26.366 

78 

U 

85.000 

86.894 

-1.894 

79-78 

-15.616 

79 

19 

70.000 

71.278 

-1.278 

80-79 

5.982 

80 

12 

72.500 

77.260 

-4.760 

81-80 

-6.458 

81 

7 

60.000 

70.802 

•10.802 

82-81 

-5.471 

10 

82 

7 

55.000 

65.331 

-10.331 

83-82 

-28.902 

11 

83 

10 

38.000 

36.429 

1.571 

84-83 

-13.586 

APPQBIX:  KEPOBT  FILES:  TIX  mtraraforMd  dita 


12         84 

13 

20.000 

22.843 

-2.843 

85-84 

7.175 

13         85 

12 

32.500 

30.018 

2.482 

86-85 

6.352 

U         86 

12 

35.000 

36.370 

■1.370 

87-86 

-1.198 

15         87 

13 

30.000 

35.172 

-5.172 

88-87 

-6.707 

16         88 

13 

30.000 

28.465 

1.535 

89-88 

.692 

17         89 

10 

33.500 

29.157 

4.343 

MEAN   TREND 

B 

52.912 

56.2S4 

MAXIMUM 

TREND  - 

110.371  a 

1973 

MINIMUM 

TREND  « 

22.843  a 

1984 

TREND  RANGE        ' 

87.528 

SPEARMANS  RHO  TREND  TEST  ON  FINAL  ANNUAL  TRENDS 

FOR  YEARS  WITH  WITH  5  OK  MORE  OBSERVATIONS 

IF  RHO  <  0  THEN  NET  TREND  IS  DOWN 

IF  RHO  >  0  THEN  NET  TREND  IS  UP 

RHO  >  -.76471   APPROX  p  *     ,000350   dof  >  15 


SSSSSSS   SEASONAL  SUMMARY  SSSSSS 


MEAN  ANNUAL  SEASONAL  MEDIAN  = 


56.548160 


SCAN  FOR  SEASONAL  EXTREMA 


LOCAL   MAX   » 

17.972840  a 

.22     » 

MAR  22 

LOCAL   MIN   * 

-11.976100  a 

.62     - 

AUG  15 

SEASONAL  MAX  '  17.972840  a  .22  >  MAR  22 
SEASONAL  MIN  *  -11.976100  8  .62  '  AUG  15 
SEASONAL  AMPLITUDE  '   29.948940 

ITERATIONS  =  5  COMPLETED 


SSS  VARIATION  IN  DE -  SEASONAL  I  ZED  RESIDUALS  SSS 

VARIATION  DUE  TO  TREND  >  70.15365  X 
VARIATION  DUE  TO  NOISE  '  29.84635  X 
APPROXIMATf  RSQ  •  70.15365  X 


MEAN  X 
STD  DEV  X 


59.029960 
30.731770 


MEAN  E 
STD  DEV  E 


1.056031 
16.824720 


ROBUST  SUMMARY  STATISTICS: 


MEDIAN  X 

56.793370 

MEDIAN  E 

.140057 

MEDIAN    iXJ 

56.793370 

MEDIAN    Ie;      < 

8.238347 

MAO  X 

23.328770 

MAC  E 

8.17U55 

.5   X   lOR 

23.284600 

.5   E    lOR 

8.595749 

R   STD  OEV 

34.587290 

R  STD  DEV       > 

12.119460 

.8.  X  ■  DE-SESEASONALIZED  RESIDUAL 

.B.  T  ■  TREND  COMPONENT 

.B.  E  >  DE-TRENDED,  DE-SESEASONALIZED  RESIDUAL 


APf>ElBIX:  REPORT  FILES:  TtX  wtramforMd  dita 


SSS     VARIATION   IN  OE-TRENDED  RESIDUALS     SSS 

VARIATION  DUE   TO  SEASONAL  =     24.71026  X 

VARIATION  DUE   TO  NOISE  >=     75.28974  X 

APPROXIMATE   RSO  x     24.71026  X 


MEAN  X 
STD  DEV  X 


53.931260 
19.232020 


MEAN  E 
STD  DEV  E 


.762146 
16.824720 


ROBUST  SUMMARY  STATISTICS: 


MEDIAN  X 
MEDIAN  iX 
MAD  X 
.5  X  lOR 
R  STD  DEV 


53.702680 
53.702680 
10.102360 
10.130100 
14.977780 


MEDIAN  E 
MEDIAN  JE| 
MAO  E 
.5  E  lOR 

R  STD  DEV 


-.153827 
8.319060 
8.17U57 
8.595750 
12.119470 


N.B.  X  E  DE-TRENDED  RESIDUAL 

N.B.  S  3  SEASONAL  COMPONENT 

N.B.  E  s  OE-TRENDED,  DE-SESEASONALIZED  RESIDUAL 


SSS  TOTAL  VARIATION  IN  SERIES:  SSS 


VARIATION  DUE  TO  TREND 

VARIATION  DUE  TO  SEASON 

VARIATION  DUE  TO  S-^T 

VARIATION  DUE  TO  NOISE 
APPROXIMATE  RSO 


58.21906  X 
-5.71735  X 
68.45300  X 
31.54700  X 
68.45300  X 


MEAN  X 
STD  DEV  X 


55.650940 
29.989430 


MEAN  E 
STD  DEV  E 


.909089 
16.824720 


ROeuST  SUMMARY  STATISTICS 


MEDIAN  X 
MEDIAN  |XJ 
MAD  X 
.5  X  lOR 
R  STD  DEV 


51.000000 
51.000000 
20.500000 
20.281250 
30.393350 


MEDIAN  E 
MEDIAN  iEi 
MAD  E 
.5  E  lOR 

R  STD  DEV 


-.006885 
8.181341 
8.17U56 
8.595751 
12.119470 


N.B.  X  «  RAW  DATA 

N.B.  T  '   TREND  COMPONENT 

N.B.  S  ■  SEASONAL  COMPONENT 

N.B.  E  -  OE-TRENDED.  DE-SESEASONALIZED  RESIDUAL 

N.B.  D  ■  CHRONOLOGICAL  SAMPLE  DATE 


SUMMARY  STATISTICS  FOR  DE-SEASONALIZED,  DE-TRENDED  RESIDUALS 
08$  212 


A -MEAN 
STD  DEV 


.9091 
16.8247 


MAX 
Q90 
075 


87.5349 
18.3060 
9.6436 


APfElOIX:  l£PORT  FILES:  TEX  ksitr«wforHd  dita 


O50 

-.0069 

025 

-7.5479 

010 

-16.4565 

MIN 

-82.4167 

RANGE  STATS 

MAX-HIN 

169.9516 

090-Q10 

34.7625 

075-025 

17.1915 

TINCEY  FENCES  (STANDARDIZED)  FOR  OUTLIER  SCAN 


UPPER 

OUTER 

FENCE  ' 

61.2181 

7.4826 

UPPER 

INNER 

FENCE  = 

35.4308 

4.3307 

LOWER 

INNER 

FENCE  ' 

-33.3352 

-4.0745 

LOWER 

OUTER 

FENCE  * 

-59.1225 

-7.2265 

Ml  OUTLIERS  IN  DE-TRENDED  t  DE- SEASONAL  I  ZED  RESIDUALS  III 
DATE  TIME  RAW  DATUM   RESIDUAL  STANDARDIZED 


RESIDUAL  730809 

1200 

162.000 

63.396 

7.749 

IS 

FAR  OUT 

RESIDUAL 

741114 

1200 

151.000 

87.535 

10.699 

IS 

FAR  OUT 

RESIDUAL 

751014 

1200 

91.000 

30.704 

3.753 

IS 

ALMOST  OUT 

RESIDUAL 

760511 

1200 

126.000 

64.465 

7.880 

IS 

FAR  OUT 

RESIDUAL 

760607 

1200 

19.000 

-35.582 

-4.349 

IS 

OUT 

RESIDUAL 

780228 

1200 

20.000 

-82.417 

-10.074 

IS 

FAR  OUT 

RESIDUAL 

790219 

1200 

121.000 

37.267 

4.555 

IS 

OUT 

RESIDUAL 

790309 

1200 

128.000 

39.319 

4.806 

IS 

OUT 

RESIDUAL 

800730 

1055 

14.000 

-51.922 

-6.346 

IS 

ALMOST  FA*  OUT 

RESIDUAL 

820914 

830 

24.000 

-30.819 

-3.767 

IS 

ALMOST  OUT 

RESIDUAL 

821018 

1130 

20.000 

-36.830 

-4.502 

IS 

OUT 

RESIDUAL 

870512 

1200 

76.000 

35.952 

4.394 

IS 

OUT 

THERE 

ARE: 

2  ALMOST  OUT  DATA 

5  OUT 

DATA 

1  ALMOST  FAR  OUT 

DATA 

4  FAR 

OUT  DATA 

12  TOTAL  OUTLIERS 

IDENTIFIED 

N.B.  OUTLIERS  ARE  STANDARDIZED  BY  DIVISION  BY  THE  MEDIAN  ABSOLUTE  RESIDUAL  >  S 

FOR  NORMAL  DISTRIBUTIONS: 
OUTER  FENCES  ARE  ABOUT  7S  OR  4.7  X  STANDARD  DEVIATION 
INNER  FENCES  ARE  ABOUT  4S  OR  2.7  X  STANDARD  DEVIATION 
ALSO  6S  OR  4  X  STANDARD  DEVIATION  (ALMOST  FAR  OUT) 
I   3.75S  OR  2.5  X  STANDARD  DEVIATION  (ALMOST  OUT) 
ARE  STANDARDIZED  VALUES  USE  BY  OTHER  AUTHORS 


TUKEY  FENCES  ARE  DEFINED  USING  1/2  THE  INTERQUARTILE  RANGE 
OF  RESIDUALS  WHICH  IS  NOT  EXACTLY  IDENTICAL  TO  S 
EXCEPT  FOR  NORMAL  DISTRIBUTIONS  BUT  IS  ALSO  A  CONMON 
ROBUST  MEASURE  OF  DISPERSION 


«PI>BBIX:  KPCRT  FILES:  TIX  log^g  trmforwd  data 


PROGRAM   TRX 

REPORT  DATA  FILE  >  O:\DILD0\TRX\CRN03L 
INPUT  DATA  FILE  -  D:\OILDO\DAT2\N03.dat 
VARIABLE  ■  N03 
{5X,i7,i5,f9.3) 
START  YEAR  ■  720000   END  YEAR  >  900000 

ZERO  OPTION  CODE  >  1 

NUMERIC  NISSING  VALUE  CODE  '         -9.00000000 

MULTIPLIER/DIVISOR  CODE  >  1 

MULTIPLIER/DIVISOR  '  .00000000 

INPUT  DATA  ARE  LEFT  AS  IS 

LOGG  CODE  =   1 

RUNS  *        1 

ITERATIONS  *       5 

MINIMUM  N  /  WINDOW  >  13 

MINIMUM  SEASONAL  SMOOTHING  WINOOU  «  1/  8  yr 

PLOT  CODES: 
PLOT  INDICATOR  '  0 
PLOT  DENSITY  1  =  0 
LOG  PLOT  TYPE  =  0 

PLOT  TYPE  CODE  =  L  (i.e.  LONG/SHORT) 
PLOT  HEADER  1  =  Nitrate  U9/L 
PLOT  HEADER  2  '   CLEARUATER  L 
RUN  DATE  HEADER  >  JUN  13  1991 
Y-AXIS  TITLE  ■  LofllO  Concentration  ug/L 
REPORT  FILE  CODE  «  SR  (i.e.  SHORT  REPORT) 
OUTPUT  FILE  CODE  '   0 
NO  OUTPUT  FILES  REQUESTED 
FIRST  LINE  OF  INPUT  DATA  FILE: 

PARAMETER  «  N03     SITE  «  CLEARUATER  L 

QUALITY  DATA  FOR  PARM  >  N03     READ  t  STORED 
212  RECORDS  WERE  READ 

THE  FOLLOWING  3  INPUT  DATA  WERE  IDENITIFIED  AS  ZERO'S: 
870901  1100      .000  ZERO 
880830  1030      .000  ZERO 
880926  1000     .000  ZERO 

SMALLEST  NON-ZERO  POSITIVE  VALUE  «     4.00000 

ZERO  OPTION  CODE  =  1 

SETTING  O'S  TO  1/2  SMALLEST  POSITIVE  NON-ZERO  VALUE  »     2.00000 


DOING  LOG10  TRANSFORM  OF  DATA  SERIES 


BEGINNING  RUN  1 


TREND  SUMMARY  ON  LOG10  SERIES 

YEAR    N    ITER  »  1  ITER  «  5 

2.006  2.084 

1.785  1.862 

1.818  1.8S3 

1.663  1.756 

1.803  1.788 

1.929  1.975 

1.845  1.855 

1.860  1.905 


1 

73 

16 

2 

74 

21 

75 

10 

76 

13 

77 

10 

78 

14 

79 

19 

8 

80 

12 

OIFF 

YR-TO-YR 

DIFF 

-.077 

74-73 

-.222 

-.077 

75-74 

-.009 

-.035 

76-75 

-.097 

-.093 

77-76 

.032 

.014 

78-77 

.186 

-.045 

79-78 

-.120 

-.010 

80-79 

.050 

-.045 

81-80 

-.039 

APPEIBIX:  KPORT  F 

IL£S: 

9         81 

7 

1.778 

.866 

10         82 

7 

1.740 

.832 

11         83 

10 

1.578 

.569 

12         8^ 

13 

1.301 

.277 

13         85 

12 

1.511 

.462 

U         86 

12 

1.5U 

.556 

15         87 

13 

1.477 

.531 

16         88 

13 

1.477 

.U9 

17         89 

10 

1.523 

.435 

MEAN  TREND   * 

1.685 

1 

.709 

MAXIMUM 

TREND 

X 

2.084  a 

1973 

MINIMUM 

TREND 

K 

1.277  S 

1984 

TREND   RANGE 

S 

.807 

»10 


tmforaad  dmtm 


.088 

82-81 

.092 

83-82 

.009 

84-83 

.024 

85-84 

.049 

86-85 

.012 

87-86 

.054 

88-87 

.028 

89-88 

.088 

-.034 
.263 
.293 
.185 
.094 
-.025 
-.082 
-.015 


TREND 

SUMMARY   ON 

NATURAL   SERIES 

YEAR 

N 

ITER   =    1 

ITER   *     5 

DIFF 

YR-TO-YR 

DIFF 

1 

73 

16 

101.499 

121.323 

-19.825 

74-73 

-48.491 

2 

74 

21 

61.000 

72.832 

-11.832 

75-74 

-1.545 

3 

75 

10 

65.810 

71.287 

-5.477 

76-75 

-14.303 

4 

76 

13 

46.000 

56.984 

-10.984 

77-76 

4.420 

5 

77 

10 

63.482 

61.403 

2.079 

78-77 

32.896 

6 

78 

14 

85.000 

94.299 

-9.299 

79-78 

-22.690 

7 

79 

19 

70.000 

71.610 

-1.610 

80-79 

8.742 

8 

80 

12 

72.457 

80.352 

-7.895 

81-80 

-6.826 

9 

81 

7 

60.000 

73.526 

-13.526 

82-81 

-5.603 

10 

82 

7 

55.000 

67.923 

-12.923 

83-82 

-30.815 

11 

83 

10 

37.881 

37.108 

.773 

84-83 

-18.187 

12 

84 

13 

20.000 

18.921 

1.079 

85-84 

10.029 

13 

85 

12 

32.404 

28.950 

3.453 

86-85 

7.016 

14 

86 

12 

35.000 

35.966 

-.966 

87-86 

-2.006 

15 

87 

13 

30.000 

33.960 

-3.960 

88-87 

-5.814 

16 

88 

13 

30.000 

28.146 

1.854 

89-88 

-.942 

17 

89 

10 

33.317 

27.204 

6.112 

G- 

MEAN  TREND   < 

48.382 

51.183 

MAXIMUM 

TREND  « 

121.323 

a  1973 

MINIMUM 

TREND   ' 

18.921 

a  1984 

TREND   RANGE        « 

102.403 

SPEARMANS  RHO  TREND  TEST  ON  FINAL  ANNUAL  TRENDS 

F0«  TEARS  WITH  UITM  5  00  MORE  OfiSERVATIONS 

IF  RHO  <  0  THEN  NET  TREND  IS  DOWN 

IF  RHO  >  0  THEN  NET  TREND  IS  UP 

RHO  «  -.77451   APPROX  p  «  .000262   dof  «  15 


SSSSSSS   SEASONAL  SUMMARY  SSSSSS 


MEAN  ANNUAL  SEASONAL  MEDIAN  ■ 
G-MEAN  ANNUAL  SEASONAL  MEDIAN  « 


1.710786  LOG  UNITS 
51.379040  NATURAL  UNITS 


SCAN  FOR  SEASONAL  EXTRENA 


LOCAL  MAX  -     .176214  a  .25  '     APR  2 
LOCAL  MIN  ■    -.106097  a   .62  >  AUG  15 


SEASONAL  MAX  = 


.176214  a  .25  =  APR  2 


APPBBIX:  KPORT  FILES:  TIX  los,.  tmforvd  (teta 


10 


SEASONAL  MIN  ■    -.106097  9   .62  >  AUG  15 
SEASONAL  AMPLITUDE  *  .282311 


IN  NATURAL  UNITS 

AT  G-MEAN  ANNUAL  SEASONAL  MEDIAN  ■ 


51.379 


SEASONAL  MAX  '  77.090340  8  .25  >  APR  2 
SEASONAL  MIN  >  40.242840  a  .62  «  AUG  15 
SEASONAL  AMPLITUDE  ■   36.847500 


AT  MAXIMUM  TREND 


«I21.323  IN  1973 


SEASONAL  MAX  =  182.036500  S  .25  <  APR  2 
SEASONAL  MIN  <  95.027050  S  .62  «  AUG  15 
SEASONAL  AMPLITUDE  '       87.009480 


AT  MINIMUM  TREND  * 


18.921   IN  1984 


SEASONAL  MAX  ^  28.389340  a  .25  >  APR  2 
SEASONAL  MIN  >  14.819860  a  .62  *  AUG  15 
SEASONAL  AMPLITUDE  -   13.569480 

ITERATIONS  >  5  COMPLETED 


SS$  VARIATION  IN  DE- SEASONAL  I ZED  RESIDUALS  $S$ 


VARIATION  DUE  TO  TREND 
VARIATION  DUE  TO  NOISE 
APPROXIMATE  RSQ 


55.40270  X 
44.59730  X 
55.40270  X 


MEAN  X 
STO  DEV  X 


1.689344 
.324584 


MEAN  E 
STD  DEV  E 


-.031989 
.214785 


ROBUST  SmWARY  STATISTICS: 


MEDIAN  X 

1.750334 

MEDIAN  E 

X 

.000181 

MEDIAN  ;x| 

1.750334 

MEDIAN  JE| 

s 

.075157 

MAO  X 

.181839 

MAD  E 

5 

.075277 

.5  X  lOR 

.196041 

.5  E  IQR 

X 

.075290 

R  STD  DEV 

.269595 

R  STD  DEV 

= 

.111606 

N.B.  X  >  DE-SESEASONALIZED  RESIDUAL 

N.B.  T  =  TREND  COMPONENT 

N.B.  E  ■  DE-TRENDED,  DE-SESEASONALIZED  RESIDUAL 


SSS  VARIATION  IN  DE-TRENDED  RESIDUALS  SSS 

VARIATION  DUE  TO  SEASONAL  «  23.71050  X 
VARIATION  DUE  TO  NOISE  e  76.28950  X 
APPROXIMATE  RSQ  m     23.71050  X 


MEAN  X 
STD  DEV  X 


1.645320 
.240116 


MEAN  E 
STD  DEV  E 


-.033648 
.214785 


APfElBIX:  KEPOBT  FILES:  TU  loo^g  tramforMd  dita 


ROBUST  SUMMARY  STATISTICS: 


MEDIAN  X 

1.664036 

MEDIAN  E 

= 

-.001478 

MEDIAN  ix| 

1.664036 

MEDIAN  |E| 

X 

.07^868 

MAO  X 

.092489 

MAD  E 

« 

.075277 

.5  X  lOR 

.092368 

.5  E  IM 

m 

.075291 

R  STD  DEV 

.137125 

R  STO  OEV 

m 

.111606 

N.B.  X  '  DE-TREMDED  RESIDUAL 

N.B.  S  >  SEASONAL  COMPONENT 

N.B.  E  '   DE-TRENDED.  DE-SESEASONALIZED  RESIDUAL 


SSS     TOTAL  VARIATION  IN  SERIES:  SSS 


VARIATION  DUE  TO 

TREND 

=  43.20377  X 

VARIATION  DUE  TO 

SEASON 

«   2.86671  X 

VARIATION  DUE  TO 

S+T 

=  56.64535  X 

VARIATION  DUE  TO 

NOISE 

=     43.35465  X 

APPROXIMATE  RSQ 

=  56.64535  X 

MEAN  X      = 

1.657526 

MEAN  E 

STO  DEV  X   = 

.325794 

STD  DEV  E 

-.032818 
.214785 


ROBUST  SUMMARY  STATISTICS 


MEDIAN  X 

1.707570 

MEDIAN  E 

-.000649 

MEDIAN  ;x; 

1.707570 

MEDIAN  |e; 

.075251 

MAO  X 

.167491 

MAD  E 

.075277 

.5  X  lOR 

.169037 

.5  E  lOR 

.075291 

R  STD  DEV 

.248323 

R  STO  DEV 

.111606 

N.B.  X  =  RAW  DATA 

N.B.  T  *   TREND  COMPONENT 

N.B.  S  '   SEASONAL  COMPONENT 

N.B.  E  s  DE-TRENOEO,  DE-SESEASONALIZEO  RESIDUAL 

N.B.  0  >  CHRONOLOGICAL  SAMPLE  DATE 


SUMPURY  STATISTICS  FOR  DE- SEASONAL  I  ZED.  DE-TRENDED  RESIDUALS 


212 


A-NEAN 

-.0328 

STD  DEV 

.2148 

NM 

.3853 

090 

.1373 

075 

.0746 

050 

-.0006 

025 

-.0760 

010 

-.1889 

MIN 

-1.1360 

RANGE  STATS 

MAX-MIN 

1.5192 

Q9O-O10 

.3263 

075-025 

.1506 

AFfODIX:  KPOT  FILES:  TIX  log^g  tmfoTMd  (teta 

TUKET  FENCES  (STANDARDIZED)  FCK  OUTLIER  SCAN 


UPPER  OUTER  FENCE  « 

UPPER  INNER  FENCE  > 

LOWER  INNER  FENCE  = 

LOWER  OUTER  FENCE  > 


.5264  6.9948 

.3005  3.9932 

-.3018  -4.0109 

-.5277  -7,0125 


HI  OUTLIERS  IN  DE- TRENDED  t  DE- SEASONAL  I  ZED  RESIDUALS  III 


DATE  TIME  RAW  DATUM 


LOG10   RESIDUAL  STANDARDIZED 


RESIDUAL 
RESIDUAL 
RESIDUAL 
RESIDUAL 
RESIDUAL 
RESIDUAL 
RESIDUAL 
RESIDUAL 
RESIDUAL 
RESIDUAL 
RESIDUAL 
RESIDUAL 
RESIDUAL 
RESIDUAL 
RESIDUAL 
RESIDUAL 
RESIDUAL 
RESIDUAL 
RESIDUAL 
RESIDUAL 
RESIDUAL 


741114 
760511 
760607 
780228 
800730 
820914 
821018 
S40306 
840828 
850924 
870512 
870805 
870901 
870928 
880803 
880830 
880926 
890516 
890613 
890925 
891128 


1200 
1200 
1200 
1200 
1055 
830 
1130 
905 
955 
950 
1200 
1000 
1100 
1010 
1020 
1030 
1000 
900 
850 
1055 
1116 


151.000 

126.000 

19.000 

20.000 

14.000 

24.000 

20.000 

50.000 

5.000 

12.000 

76.000 

13.000 

2.000 

4.000 

5.000 

2.000 

2.000 

57.000 

51.000 

4.OO0 

10.000 


2.179 

2.100 

1.279 

1.301 

1.146 

1.380 

1.301 

1.699 

.699 

1.079 

1.881 

1.114 

.301 

.602 

.699 

.301 

.301 

1.756 

1.708 

.602 

1.000 


.383 

.313 

-.446 

-.788 

-.660 

-.366 

-.452 

.300 

-.481 

-.300 

.319 

-.315 

-1.136 

-.848 

-.649 

■1.053 

■1.067 

.302 

.314 

-.751 

-.383 


5.093 

4.155 

-5.932 

-10.477 

-8.768 

-4.862 

-6.011 

3.981 

-6.391 

-3.989 

4.240 

-4.180 

-15.096 

■11.263 

-8.622 

■13.993 

•14.175 

4.007 

4.167 

-9.974 

-5.090 


IS  OUT 

IS  OUT 

IS  OUT 

IS  FAR  OUT 

IS  FAR  OUT 

IS  OUT 

IS  ALMOST  FAR  OUT 

IS  ALMOST  OUT 

IS  ALMOST  FAR  OUT 

IS  ALMOST  OUT 

IS  OUT 

IS  OUT 

IS  FAR  OUT 

IS  FAR  OUT 

IS  FAR  OUT 

IS  FAR  OUT 

IS  FAR  OUT 

IS  OUT 

IS  OUT 

IS  FAR  OUT 

IS  OUT 


THERE  ARE: 


2  ALMOST- OUT  DATA 
9  OUT  DATA 

2  ALMOST  FAR  OUT  DATA 
8  FAR  OUT  DATA 
21  TOTAL  OUTLIERS  IDENTIFIED 


N.B.  OUTLIERS  ARE  STAWARDIZEO  BY  DIVISION  BY  THE  >£DIAN  ABSOLUTE  RESIDUAL  =  S 

FOR  NORMAL  DISTRIBUTIONS: 
OUTER  FENCES  ARE  ABOUT  7S  OR  4.7  X  STANDARD  DEVIATION 
INNER  FENCES  ARE  ABOUT  4S  OR  2.7  X  STANDARD  DEVIATION 
ALSO  6S  OR  4  X  STANDARD  DEVIATION  (ALMOST  FAR  OUT) 
t  3.75S  OR  2.5  X  STANDARD  DEVIATION  (ALMOST  OUT) 
ARE  STANDARDIZED  VALUES  USE  BY  OTHER  AUTHORS 

N.B.  TUKET  FENCES  ARE  DEFINED  USING  1/2  THE  INTERQUARTILE  RANGE 
OF  RESIDUALS  WHICH  IS  NOT  EXACTLY  IDENTICAL  TO  S 
EXCEPT  FOK  NORMAL  DISTRIBUTIONS  BUT  IS  ALSO  A  CONMON 
ROBUST  MEASURE  OF  DISPERSION 


WPODIX:  lEPaiT  FILES:  Til  iwttnrmformmd  dit«  10 


PROGRAM  TR1 

REPORT  DATA  FILE  *  D:\DlL00\TR1\CRno3 
INPUT  DATA  FILE  >  D:\DILDO\DAT2\N03.dat 
VARIABLE  *   N03 
(5X,)7,i5,f9.3) 
START  YEAR  ■  720000   END  YEAR  ■  920000 

ZERO  OPTION  CODE  *   1 

NUMERIC  MISSING  VALUE  CODE  >    -9.00000000 

MULTIPLIER/DIVISOR  CODE  =  1 

MULTIPLIER/DIVISOR  '  .00000000 

INPUT  DATA  ARE  LEFT  AS  IS 

LOGG  CODE  =       0 

RUNS  *        2 

FILTER  =  5 

SAMPLING  INTERVAL  CODE  '  1 

MINIMUM  SMOOTHING  BASED  ON  MEDIAN  SAMPLING  INTERVAL 

UNIT  SMOOTHING  INTERVAL  •   1.00  yr 

PLOT  CWES: 
PLOT  INDICATOR  =   0 
PLOT  DENSITY  =  D 
LOG  PLOT  TYPE  =  0 

PLOT  TYPE  CODE  *  I     (i.e.  LONG/SHORT) 
PLOT  HEADER  1  "   Nitrate  ug/L 
PLOT  HEADER  2  =  CLEARWATER  L 
RUN  DATE  HEADER  =   JUN  13  1991 
Y-AXIS  TITLE  «  Concentration  ug/L 
REPORT  FILE  CODE  *   SR  (i.e.  SHORT  REPORT) 
OUTPUT  FILE  CODE  =  0 
NO  OUTPUT  FILES  REQUESTED 
FIRST  LINE  OF  INPUT  DATA  FILE: 

PARAMETER  *   N03      SITE  *   CLEARUATER  L 

QUALITY  DATA  FOR  PARM  <  N03     READ  t  STORED 
213  RECORDS  WERE  READ 

THE  FOLLOWING  3  INPUT  DATA  WERE  IDENITIFIEO  AS  ZERO'S: 
870901  1100      .000  ZERO 
880830  1030      .000  ZERO 
880926  1000      .000  ZERO 

SMALLEST  NON-ZERO  POSITIVE  VALUE  >     4.00000 

ZERO  OPTION  CODE  *   1 

SETTING  O'S  TO  1/2  SMALLEST  POSITIVE  NON-ZERO  VALUE  ■     2.00000 


BEGINNING  RUN  1 

UNIT  SMOOTHING  INTERVAL  (USD  =      1.00000  yr 

LENGTH  OF  RECORD  e  16.60808  yr 

MEAN  SAMPLING  INTERVAL  >       .07834  yr 

MEAN  SAMPLING  RATE/US  1  =  12.76487  #/USI 
MEDIAN  SAMPLING  INTERVAL  =       .07382  yr 

MEDIAN  SAMPLING  RATE/USI  *  13.S4609  «/USI 

MINIMUM  SMOOTHING  '  14  OBSERVATIONS 


TREND  SUMMARY: 
GRAND  TIME -WEIGHTED  MEAN  TREND  «     S3.3S4 
SCAN  SMOOTH  FOR  LOCAL  EXTREMA  INCLUDING  END  POINTS 


APf>EIDIX:  KPCRT  FILES:  Til  intrwrmformtd  dita 


LOCAL  MAX  < 
LOCAL  MIN  < 
LOCAL  MAX  < 
LOCAL  MIN  = 
LOCAL  MAX  * 
LOCAL  MIN  = 
LOCAL  MAX  ' 
LOCAL  MIN  > 
LOCAL  MAX  ■ 
LOCAL  MIN  : 
LOCAL  MAX  = 
LOCAL  MIN  = 
LOCAL  MAX  r 
LOCAL  MIN  < 

MAX  TREND  ' 
MIN  TREND  ' 
TREND  RANGE 


110.130 
66.^06 
67.009 
55.109 
87.733 
74.436 
75.176 
24.383 
31.637 
29.892 
39.577 
26.131 
28.713 
10.904 


73.489 
74.788 
74.870 
76.818 
78.656 
79.793 
80.362 
84.657 
85.387 
85.807 
87.113 
88.736 
89.223 
90.097 


JUN  28 
OCT  15 
NOV  14 
OCT  26 
AUG  28 
OCT  17 
MAY  12 
AUG  28 
MAY  22 
OCT  22 
FEB  11 
SEP  26 
MAR  23 
FEB  5 


1973 
1974 
1974 
1976 
1978 
1979 
1980 
1984 
1985 
1985 
1987 
1988 
1989 
1990 


110.130  a  73.489  »  JUN  28  1973 
10.904  S  90.097  «  FEB  5  1990 
99.225 


ANNUAL   TIME-WEIGHTED  MEAN   TREND   SUMMARY 

YEAR 

N 

MEAN 

TR-TO-YR 

OIFF 

1 

73 

16 

104.811 

74-73 

-36.470 

2 

74 

21 

68.342 

75-74 

-2.597 

3 

75 

10 

65.744 

76-75 

-8.198 

4 

76 

13 

57.546 

77-76 

8.456 

5 

77 

10 

66.003 

78-77 

18.847 

6 

78 

14 

84.850 

79-78 

-7.236 

7 

79 

19 

77.613 

80-79 

-3.276 

8 

80 

12 

74.338 

81-80 

-7.771 

9 

81 

7 

66.567 

82-81 

-13.120 

10 

82 

7 

53.U7 

83-82 

-16.240 

11 

83 

10 

37.207 

84-83 

-10.984 

12 

84 

13 

26.223 

85-84 

4.390 

13 

85 

12 

30.614 

86-85 

4.354 

14 

86 

12 

34.968 

87-86 

-.074 

15 

87 

13 

34.894 

88-87 

-7.294 

16 

88 

13 

27.600 

89-88 

-3.480 

17 

89 

10 

24.120 

90-89 

-13.216 

18 

90 

1 

10.904 

MEAN 

1  TREND 

= 

52. 5U 

MAXIMUM  TREND   « 

104.811  a 

1973 

MINIMUM   TREND   > 

10.904  a 

1990 

TREND  RANGE        > 

93.907 

SPEARMANS  RHO  TREND  TEST  (2-SIDED)  ON  FINAL  ANNUAL  TRENDS 

FOR  YEARS  WITH  WITH  5  OR  MOKE  OBSERVATIONS 

IF  RHO  <  0  THEN  MET  TREND  IS  DOWN 

IF  RHO  >  0  THEN  MET  TREND  IS  UP 

RHO  s  -.82108   APPROXX  p  <  .000054   dof  -     15 


SSS  APPROXIMATE  SOURCES  OF  VARIATION  SSS 


VARIATION  DUE   TO  TREND 

62.67897  X 

VARIATION  DUE   TO  NOISE 

37.32103  X 

APPROXIMATE   RSO 

62.67897  X 

MEAN  X             =         55.469480 

MEAN  E             = 

-1.075196 

APPQBIX:  KEPORT  FIL£S:  Til  tMttrt^miorwad  dmtM 
30.035590  STO  OEV  E       =         18.363080 
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ROBUST  SUMMARY  STATISTICS: 


MEDIAN  X    « 

51.000000 

MEDIAN  E 

-1.289604 

MEDIAN  ix|   « 

51.000000 

MEDIAN  |E| 

9.891596 

MAD  X       ■ 

21.000000 

MAD  E 

9.116722 

.5  X  lOR 

20.812500 

.5  E  lOR 

9.372667 

R  STD  DEV   « 

31.134650 

R  STD  DEV 

13.516470 

N.B.  X  >  RAW  DATA 

N.B.  T  <  TREND 

COMPONENT 

N.B.  E  X  RESIDUAL  =  X-T 

SUWURY  STATISTICS  FOR  RESIDUALS 
oes  213 


A-HEAN 

-1.0752 

STD  DEV 

18.3631 

MAX 

83.9915 

090 

19.8664 

075 

7.0706 

050 

-1.2896 

025 

-11.6748 

010 

-22.7806 

NIN 

-61.2132 

RANGE  STATS 

MAX-HIN 

145.2046 

090-010 

42.6469 

075-025 

18.7453 

TUKEY  FENCES  (STANDARDIZED)  FOR  OUTLIER  SCAN 


UPPER 

OUTER 

FENCE  - 

63.3066 

6.4000 

UPPER 

INNER 

FENCE  = 

35.1886 

3.5574 

LOWER 

INNER 

FENCE  * 

-39.7928 

-4.0229 

LOWER 

OUTER 

FENCE  - 

-67.9108 

-6.8655 

III  OUTLIERS  IN  RESIDUALS  III 


DATE  TIME  RAW  DATUM   RESIDUAL  STANDARDIZED 


RESIDUAL 
RESIDUAL 
RESIDUAL 
RESIDUAL 
RESIDUAL 
RESIDUAL 
RESIDUAL 
RESIDUAL 
RESIDUAL 
RESIDUAL 
RESIDUAL 


730809  1200 
741114  1200 
760429  1200 
760511  1200 
760607  1200 
780228  1200 
790219  1200 
790309  1200 
790710  930 
800730  1055 
870512  1200 


162.000 

151.000 

96.000 

126.000 

19.000 

20.000 

121.000 

128.000 

35.000 

14.000 

76.000 


57.019 
83.991 
36.050 
66.579 
-39.230 
-61.213 
38.247 
46.090 
-40.745 
-60.295 
38.788 


1  ALMOST  OUT  DATA 
6  OUT  DATA 

2  ALMOST  FAR  CUT  DATA 


5.764 
8.491 
3.6U 
6.731 
-3.966 
•6.188 
3.867 
4.659 
-4.119 
-6.096 
3.921 


IS  OUT 

IS  FAR  OUT 

IS  OUT 

IS  FAR  OUT 

IS  ALMOST  OUT 

IS  ALMOST  FAR  OUT 

IS  OUT 

IS  OUT 

IS  OUT 

IS  ALMOST  FAR  OUT 

IS  OUT 


APPQDIX:  KEPOtT  FILES:  Til  wtrwwfoTMd  data  13 

2  FAR  CUT  DATA 
11  TOTAL  OUTLIERS  IDENTIFIED 


N.B.  OUTLIERS  ARE  STANDARDIZED  BY  DIVISION  ST  THE  MEDIAN  ABSOLUTE  RESIDUAL  r  S 

FOR  NORMAL  DISTRIBUTIONS: 
OUTER  FENCES  ARE  ABOUT  7S  OR  4.7  X  STANDARD  DEVIATION 
INNER  FENCES  ARE  ABOUT  4S  OR  2.7  X  STANDARD  DEVIATION 
ALSO  65  OR  4  X  STANDARD  DEVIATION  (ALMOST  FAR  OUT) 
t  3.75S  OR  2.5  X  STANDARD  DEVIATION  (ALMOST  OUT) 
ARE  STANDARDIZED  VALUES  USE  BT  OTHER  AUTHORS 

N.B.   TUKEY  FENCES  ARE  DEFINED  USING  1/2  THE  INTERQUARTILE  RANGE 
OF  RESIDUALS  WHICH  IS  HOT  EXACTLY  IDENTICAL  TO  S 
EXCEPT  FOR  NORMAL  DISTRIBUTIONS  BUT  IS  ALSO  A  OM«0N 
ROBUST  MEASURE  OF  DISPERSION 


M>PEieiX:  lEPOIT  FILES:  Til   los^g  tnraforMd  dita 


PROGRAM  TR1 

REPORT  DATA  FILE  ■  D:\DILD0\TR1\CRno3L 
INPUT  DATA  FILE  >  D:\0JL0O\0AT2\NO3. dat 
VARIABLE  ■  N03 
{5X,i7.15.f9.3) 
START  YEAR  »  720000   END  YEAR  «  920000 

ZERO  OPTION  CODE  *   1 

NLMERIC  MISSING  VALUE  CODE  '    -9.00000000 

MULTIPLIER/DIVISOR  CODE  ■  1 

MULTIPLIER/DIVISOR  =      .00000000 

INPUT  DATA  ARE  LEFT  AS  IS 

LOCG  CODE  >   1 

RUNS  '   1 

FILTER  «  5 

SAMPLING  INTERVAL  CODE  ■  1 

MINIMUM  SMOOTHING  BASED  ON  MEDIAN  SAMPLING  INTERVAL 

UNIT  SMOOTHING  INTERVAL  =   1.00  yr 

PLOT  COOES: 
PLOT  INDICATOR  *   0 
PLOT  DENSITY  =  D 
LOG  PLOT  TYPE  =  0 

PLOT  TYPE  CODE  =  L   (i.e.  LONG/SHORT) 
PLOT  HEADER  1  <  Nitrate  ug/L 
PLOT  HEADER  2  >  CLEARUATER  L 
RUN  DATE  HEADER  '   JUN  13  1991 
Y-AXIS  TITLE  *   LoglO  Concentration  ug/L 
REPORT  FILE  CODE  >  SR  (i.e.  SHORT  REPORT) 
OUTPUT  FILE  CODE  «  0 
NO  OUTPUT  FILES  REQUESTED 
FIRST  LINE  OF  INPUT  DATA  FILE: 

PARAMETER  <  N03      SITE  ■  CLEARWATER  L 

QUALITY  DATA  FOR  PARM  >  N03     READ  I   STORED 
213  RECORDS  WERE  READ 

THE  FOLLOWING  3  INPUT  DATA  WERE  IDENITIFIED  AS  ZERO'S: 
870901  1100      .000  ZERO 
880830  1030     .000  ZERO 
880926  1000      .000  ZERO 

SMALLEST  NON-ZERO  POSITIVE  VALUE  >     4.00000 

ZERO  OPTION  CODE  *   1 

SETTING  O'S  TO  1/2  SMALLEST  POSITIVE  NON-ZERO  VALUE  >     2.00000 


DOING  LOG10  TRANSFORM  OF  DATA  SERIES 


BEGINNING  RUN  1 


UNIT  SMOOTHING  INTERVAL  (USD  '  1.00000  yr 


LENGTH  OF  RECORD 
MEAN  SAMPLING  INTERVAL 
MEAN  SAMPLING  RATE/US  I 
MEDIAN  SAMPLING  INTERVAL 
MEDIAN  SAMPLING  RATE/USl 
MINIMUM  SMOOTHING 


TREND  SUMMARY: 


16.60808  yr 
.07834  yr 

12.76487  #/USI 
.07382  yr 

13.54609  #/USI 
14  OBSERVATIONS 


APfiQSIX:  lEPORT  FILES:  Ttl   log^g  trmforMd  dit«  15 

GRAND  TIME-WEIGHTED  MEAN  TREND  '  1.676     LOG  UNITS 

47.409     NATURAL  UNITS 

SCAN  SMOOTH   FOR   LOCAL   EXTREMA   INCLUDING  END  POINTS 

LOCAL  MAX  s  113.866  3  73.469  '  JUN  28  1973 

LOCAL  MIN  *  66.730  8  74.368  ■  MAY   15  1974 

LOCAL  MAX  >  71.272  8  74.870  «  NOV  14  1974 

LOCAL  MIN  -  S4.671  8  76.818  >  OCT  26  1976 

LOCAL  MAX  >  88.183  8  78.6S6  >  AUG  28  1978 

LOCAL  MIN  *  71.584  8  79.637  ■  AUG  21  1979 

LOCAL  MAX  *  73.838  8  80.362  ■  MAT  12  1980 

LOCAL  MIN  -  22.910  8  84.507  ■  JUL     4  1984 

LOCAL  MAX  :  29.240  8  85.483  «  JUN  26  1985 

LOCAL  MIN  :  28.836  8  85.807  «  OCT  22  1985 

LOCAL  MAX  s  38.663  8  87.113  =  FEB  11  1987 

LOCAL  MIN  c  11.456  8  90.097  ■  FEB     5  1990 

MAX  TREW>  >         113.886  8     73.489  ■  JUN  28  1973 
MIN  TREND  s  11.456  8     90.097  >  FEB     5   1990 

TREND  RANGE   *  102.430 


ANNUAL   TIME-WEIGHTED  MEAN   TREND  SUMMARY  ON   LOG10  SER    lES 


73 

16 

74 

21 

75 

10 

76 

13 

77 

10 

78 

14 

79 

19 

80 

12 

81 

7 

82 

7 

83 

10 

84 

13 

85 

12 

86 

12 

87 

13 

88 

13 

89 

10 

90 

1 

MEAN 

TREND 

MAXIMUM 

TREND 

MINIMUM 

TREND 

TREND  RANGE 

MEAN 

YR-TO-YR 

OIFF 

2.027 

74-73 

-.194 

1.833 

75-74 

-.004 

1.828 

76-75 

-.065 

1.764 

77-76 

.046 

1.809 

78-77 

.119 

1.928 

79-78 

-.050 

1.878 

80-79 

-.014 

1.864 

81-80 

-.042 

1.822 

82-81 

-.088 

1.734 

83-82 

-.178 

1.557 

84-83 

-.172 

1.385 

85-84 

.072 

1.457 

86-85 

.072 

1.528 

87-86 

.015 

1.543 

88-87 

-.096 

1.447 

89-88 

-.141 

1.306 

90-89 

-.247 

1.059 

1.654 

2.027  8  1973 
1.059  8  1990 
.968 


ANNUAL  T-W  MEAN  TREND  SUMMARY  ON  NATURAL  SERIES 


YEAR 

N 

MEAN 

YR-TO-YR 

DIFF 

1 

73 

16 

106.395 

74-73 

-38.366 

2 

74 

21 

68.029 

75-74 

-.694 

3 

75 

10 

67.335 

76-75 

-9.310 

4 

76 

13 

58.025 

77-76 

6.427 

5 

77 

10 

64.452 

78-77 

20.246 

6 

78 

14 

84.698 

79-78 

-9.161 

7 

79 

19 

75.537 

80-79 

-2.368 

8 

80 

12 

73.169 

81-80 

-6.767 

9 

81 

7 

66.402 

82-81 

-12.142 

10 

82 

7 

54.260 

83-82 

-18.238 

APPEIDIX:  KEPCRT  FIUS:  Til    log^g  t"^"*****^  «*•*■ 
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11    83 

10 

36.022 

84-83 

-11.775 

12    8A 

13 

24.247 

85-84 

4.382 

13    85 

12 

28.629 

86-85 

5.124 

U    86 

12 

33.754 

87-86 

1.172 

15    87 

13 

34.926 

88-87 

-6.904 

16    88 

13 

28.022 

89-88 

-7.783 

17    89 

10 

20.239 

90-89 

-8.783 

18    90 

1 

20.239 

G-MEAN  TREND  - 

47.409 

MAXIMUM 

TREND  ■ 

106.395  a  1973 

MINIMUM 

TREND  « 

11.456  8 

1990 

TREND  RANGE   ' 

94.939 

SPEARMANS  RHO  TREND  TEST  (2-SIDED)  ON  FINAL  ANNUAL  TRENDS 

FOR  YEARS  WITH  WITH  5  OR  MORE  OBSERVATIONS 

IF  RHO  <  0  THEN  MET  TREND  IS  DOUN 

IF  RHO  >  0  THEN  MET  TREND  IS  UP 

RHO  «  -.85824   APPROXX  p  =   .000026   dot  *      15 


SSS  APPROXIMATE  SOURCES  OF  VARIATION  S%% 


VARIATION  DUE  TO  TREND 
VARIATION  DUE  TO  NOISE 
APPROXIMATE  RSO 


48.30162  X 
51.69838  X 
48.30162  X 


MEAN  X 
STD  DEV  X 


1.655521 
.326339 


MEAN  E 
STD  DEV  E 


.045038 
.230726 


ROBUST  SUMMARY  STATISTICS: 


MEDIAN  X 
MEDIAN  |X| 
MAD  X 
.5  X  IQR 
R  STD  DEV 


1.707570 

1.707570 

.167491 

.175905 

.248323 


N.B.  X  ■  RAW  DATA 

N.B.  T  »  TREND  COMPONENT 

N.B.  E  ■  RESIDUAL  «  X-T 


MEDIAN  E 
MEDIAN  \l\ 
MAD  E 
.5  E  IQR 
R  STD  DEV 


.007329 
.086720 
.087448 
.095139 
.129650 


SUMMARY  STATISTICS  FOR  RESIDUALS 
OBS  213 


A-MEAN 

-.0450 

STD  DEV 

.2307 

MAX 

.3836 

Q90 

.1712 

Q75 

.0718 

Q50 

-.0073 

Q25 

-.1185 

Q10 

-.2303 

MIH 

-1.2188 

RANGE  STATS 

MAX-MIN 

1.6023 

UPKmiX:  KPQIT  FILES:  Ttl   log^g  trwrnforami  (teta 
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090-010 
075-025 


.4015 
.1903 


TUKEY  FENCES  (STANDARDIZED)  FOR  OUTLIER  SCAN 


UPPER 

OUTER 

FENCE  ' 

.6426 

7.4104 

UPPER 

INNER 

FENCE  = 

.3572 

4.1192 

LOUER 

INNER 

FENCE  * 

- .4039 

-4.6574 

LOWER 

OUTER 

FENCE  ■ 

-.6893 

-7.9486 

III  OUTLIERS  IN  RESIDUALS  III 


DATE  TINE  RAW  DATUM 


LOG10   RESIDUAL  STANDARDIZED 


RESIDUAL 

741114 

1200    151.000 

2.179 

.326 

3.760 

IS 

ALMOST  OUT 

RESIDUAL 

760607 

1200     19.000 

1.279 

-.494 

-5.691 

IS 

OUT 

RESIDUAL 

780??8 

1200    20.000 

1.301 

-.604 

-6.971 

IS 

ALMOST  FAR  OUT 

RESIDUAL 

800730 

1055     14.000 

1.146 

-.719 

-8.293 

IS 

FAR  OUT 

RESIDUAL 

820914 

830     24.000 

1.380 

-.334 

-3.853 

IS 

ALMOST  OUT 

RESIDUAL 

821018 

1130    20.000 

1.301 

-.396 

-4.570 

IS 

ALMOST  OUT 

RESIDUAL 

840828 

955     5.000 

.699 

-.662 

-7.635 

IS 

ALMOST  FAR  OUT 

RESIDUAL 

840925 

920     10.000 

1.000 

-.371 

-4.283 

IS 

ALMOST  OUT 

RESIDUAL 

850924 

950     12.000 

1.079 

-.382 

-4.407 

IS 

ALMOST  OUT 

RESIDUAL 

870805 

1000     13.000 

1.114 

-.418 

-4.820 

IS 

OUT 

RESIDUAL 

870901 

1100     2.000 

.301 

-1.219 

-14.054 

IS 

FAR  OUT 

RESIDUAL 

870928 

1010     4.000 

.602 

-.908 

-10.467 

IS 

FAR  OUT 

RESIDUAL 

880706 

1010     13.000 

1.114 

-.330 

-3.805 

IS 

ALMOST  OUT 

RESIDUAL 

R«n«n3 

1020     5.000 

.699 

-.740 

-8.528 

IS 

FAR  OUT 

RESIDUAL 

880830 

1030     2.000 

.301 

-1.132 

-13.056 

IS 

FAR  OUT 

RESIDUAL 

880926 

1000      2.000 

.501 

-1.128 

-13.002 

IS 

FAR  OUT 

RESIDUAL 

890516 

900     57.000 

1.756 

.384 

4.423 

IS 

OUT 

RESIDUAL 

890613 

850     51.000 

1.708 

.360 

4.147 

IS 

OUT 

RESIDUAL 

890925 

1055     4.000 

.602 

-.627 

-7.225 

IS 

ALMOST  FAR  OUT 

THERE 

ARE: 

6  ALMOST  OUT  DATA 

4  OUT  DATA 

3  ALMOST  FAR  OUT 

DATA 

6  FAR  OUT  DATA 

19  TOTAL  OUTLIERS 

IDENTIFIED 

N.B.  (UTLIPRS  ARE  STANDARDIZED  BY  DIVISION  BY  THE  MEDIAN  ABSOLUTE  RESIDUAL  >  S 

FOR  NORMAL  DISTRIBUTIONS: 
OUTER  FENCES  ARE  ABOUT  7S  OR  4.7  X  STANDARD  DEVIATION 
INNER  FENCES  ARE  ABOUT  4S  OR  2.7  X  STANDARD  DEVIATION 
ALSO  6S  OR  4  X  STANDARD  DEVIATION  (ALMOST  FAR  OUT) 
t   3.75S  OR  2.5  X  STANDARD  DEVIATION  (ALMOST  OUT) 
ARE  STANDARDIZED  VALUES  USE  BY  OTHER  AUTHORS 

N.B.  TUKEY  FENCES  ARE  DEFINED  USING  1/2  THE  INTEROUARTILE  MUX 
OF  RESIDUALS  WHICH  IS  NOT  EXACTLY  IDENTICAL  TO  S 
EXCEPT  FOR  NORMAL  DISTRIBUTIONS  BUT  IS  ALSO  A  COMMON 
ROBUST  MEASURE  OF  DISPERSION 


APPQBIX:  tEPORT  FILES:  TU  witmfoTHd  dita  18 


PROGRAM  TR2 

REPORT  DATA  FILE  ■  O:\DILDO\TR2\CRno3 
INPUT  DATA  FILE  «  D:\DIL0O\DAT2\N03.dat 
VARIABLE  *   1(03 
(5X,i7,J5,f9.3) 
START  YEAR  <  720000   EMO  TEAR  ■  920000 

ZERO  OPTION  CODE  >  1 

WJHERIC  MISSING  VALUE  CODE  '         -9.00000000 

NULTIPLIER/DIVISOR  CCOE  <  1 

MULTIPLIER/DIVISOR  *  .00000000 

INPUT  DATA  ARE  LEFT  AS  IS 

LOGG  CODE  >   0 

RUNS  <   1 

ITERATIONS  '       5 

FILTER  -  5 

SAMPLING  INTERVAL  CODE  *      1 

MINIMUM  SMOOTHING  BASED  ON  MEDIAN  SAMPLING  INTERVAL 

UNIT  SMOOTHING  INTERVAL  ■   1.00  yr 

MINIMUM  N  /  WINDOW  *      13 

MINIMUM  SEASONAL  SMOOTHING  WINDOW  =1/  8  yr 

PLOT  COOES: 
PLOT  INDICATOR  =  0 
PLOT  DENSITY  1  =  D 
LOG  PLOT  TYPE  *   0 

PLOT  TYPE  CODE  =  L  (i.e.  LONG/SHORT) 
PLOT  HEADER  1  '   Nitrate  ug/L 
PLOT  HEADER  2  >  CLEARWATER  L 
RUN  DATE  HEADER  >  JUN  13  1991 
Y-AXIS  TITLE  >  Concentration  ug/L 
REPORT  FILE  CODE  «  SR  (i.e.  SHORT  REPORT) 
OUTPUT  FILE  CODE  «  0 
NO  OUTPUT  FILES  REQUESTED 
FIRST  LINE  OF  INPUT  DATA  FILE: 

PARAMETER  '   N03      SITE  «  CLEARWATER  L 

QUALITY  DATA  FOR  PARM  e  n03     READ  t  STORED 
213  RECORDS  WERE  READ 

THE  FOLLOWING  3  INPUT  DATA  WERE  IDENITIFIEO  AS  ZERO'S: 
870901  1100      .000  ZERO 
880830  1030     .000  ZERO 
880926  1000      .000  ZERO 

SMALLEST  NON-ZERO  POSITIVE  VALUE  >     4.00000 

ZERO  OPTION  CtSE  >  1 

SETTING  O'S  TO  1/2  SMALLEST  POSITIVE  NON-ZERO  VALUE  '     2.00000 


BEGINNING  RUN   1 

UNIT  SMOOTHING  INTERVAL  (USD  '  1.00000  yr 

LENGTH  OF  RECORD        «  16.60808  yr 
MEAN  SAMPLING  INTERVAL  *  .07834  yr 

MEAN  SAMPLING  RATE/USI   ■  12.76487  «/USI 
MEDIAN  SAMPLING  INTERVAL  ■       .07382  yr 

MEDIAN  SAMPLING  RATE/USI  -  13.54609  «/USI 

MINIMUM  SMOOTHING        >  14  OBSERVATIONS 


TREW)  SUWtARY: 


APPQBIX:  lEPORT  FILES:  Tt2  witnnsforMd  dita 
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GRAMD  T  IK -WEIGHTED  MEAN  TREND  < 


55.830 


SCAN  SMOOTH  FOR  LOCAL  EXTREMA  INCLUDING  END  POINTS 


LOCAL  MAX  r 
LOCAL  MIN  s 
LOCAL  MAX  * 
LOCAL  MIN  « 
LOCAL  MAX  « 
LOCAL  MIN  ■ 
LOCAL  MAX  * 
LOCAL  MIN  * 
LOCAL  MAX  : 
LOCAL  MIN  * 
LOCAL  MAX  ' 
LOCAL  MIN  * 
LOCAL  MAX  < 
LOCAL  MIN  * 

MAX  TREND  ■ 
MIN  TREND  * 
TREND  RANGE 


122.408  a 
74.368  a 
74.900  8 
58.842  a 
91.517  a 
76.350  a 
76.979  a 
25.084  a 
32.026  a 
30.215  a 
38.758  a 
25.412  a 
27.214  a 
11.397  a 


73.489 
74.788 
74.870 
76.818 
78.656 
79.793 
80.362 
84.507 
85.387 
85.946 
86.960 
88.736 
89.223 
90.097 


122.408  a  73.489 
11.397  a  90.097 
111.011 


JUN  28 
OCT  15 
NOV  14 
OCT  26 
AUG  28 
OCT  17 
MAY  12 
JUL  4 
MAT  22 
DEC  12 
DEC  17 
SEP  26 
MAR  23 
FEB  5 


1973 
1974 
1974 
1976 
1978 
1979 
1980 
1984 
1985 
1985 
1986 
1988 
1989 
1990 


JUN  28  1973 
FEB  5  1990 


ANNUAL   TIME-WEIGHTED  MEAN 

TREND  SUMMARY 

YEAR 

N          ITER  *   1 

ITER   *     5 

DIFF 

YR-TO-YR 

DIFF 

1 

73 

16           104.811 

115.881 

-11.070 

74-73 

-39.032 

2 

74 

21             68.342 

76.849 

-8.507 

75-74 

-6.731 

3 

75 

10            65. 7U 

70.118 

-4.373 

76-75 

-7.964 

4 

76 

13            57.546 

62.154 

-4.607 

77-76 

6.461 

5 

77 

10            66.003 

68.615 

-2.612 

78-77 

19.653 

6 

78 

14            84.850 

88.268 

-3.418 

79-78 

-8.404 

7 

79 

19            77.613 

79.864 

-2.251 

80-79 

-3.115 

8 

80 

12            74.338 

76.750 

-2.412 

81-80 

-4.585 

9 

81 

7            66.567 

72.165 

-5.598 

82-81 

-14.698 

10 

82 

7             53. U7 

57.466 

-4.019 

83-82 

-19.504 

11 

83 

10            37.207 

37.962 

-.755 

84-83 

-11.384 

12 

84 

13             26.223 

26.578 

-.355 

85-84 

4.423 

13 

85 

12            30.614 

31.001 

-.388 

86-85 

3.902 

14 

86 

12            34.968 

34.903 

.065 

87-86 

-1.074 

15 

87 

13            34.894 

33.829 

1.064 

88-87 

-7.503 

16 

88 

13            27.600 

26.326 

1.274 

89-88 

-3.165 

17 

89 

10            24.120 

23.161 

.959 

90-89 

-11.764 

18 

90 

1             10.904 

11.397 

-.493 

MEAN 

TREND 

53.354 

55.183 

MAXIMUM  TREND  «  115.881  a  1973 
MINIMUM  TREND  «  11.397  a  1990 
TREND  RANGE   ■   104.484 

SPEARMAHS  RHO  TREND  TEST  ON  FINAL  ANNUAL  TRENDS 

FOR  YEARS  WITH  WITH  5  OR  MORE  OBSERVATIONS 

IF  RHO  <  0  THEN  MET  TREND  IS  DOWN 

IF  RHO  >  0  THEN  MET  TREND  IS  [» 

RHO  =  -.85049   APPROXX  p  »  .000015   dof  =  15 


SSSSSU   SEASONAL  SUMMARY  SSSSSS 
MEAN  ANNUAL  SEASONAL  MEDIAN  > 


S5. 093070 


U>VODlti  KPORT  FILES:  TI2  mtrmforwd  data 
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SCAN  FOit  SEASONAL  EXTREKA 


LOCAL  MAX 
LOCAL  MIN 


17.941220  a 
-13.935940  a 


.21 
.62 


SEASONAL  NAX  =  17.941220  a  .21 
SEASONAL  MIN  ^  -13.935940  8  .62 
SEASONAL  AHPLITUOE  -   31.877160 


NAR  18 
AUG  15 


NAR  18 
AUG  15 


ITERATIONS  =  5  COMPLETED 


SSS  VARIATION  IN  DE- SEASONAL  I  ZED  RESIDUALS  tSS 


VARIATION  DUE  TO  TREND 
VARIATION  DUE  TO  NOISE 
APPROXIMATE  RSO 


73.10892  X 
26.89107  X 
73.10892  X 


MEAN  X 
STO  DEV  X 


59.081190 
30.953310 


MEAN  E 
STD  DEV  E 


-.426665 
16.134380 


ROeUST  SUMMARY  STATISTICS: 


MEDIAN  X 
MEDIAN  |X 
MAD  X 
.5  X  lOR 
R  STO  DEV 


55.701320 
55.701320 
23.790020 
23.812310 
35.271130 


MEDIAN  E 

MEDIAN  JE| 
MAD  E 
.5  E  lOR 
R  STD  DEV 


-.555134 
7.743704 
7.283405 
7.597422 
10.798390 


N.B.  X  *   DE-SE SEASONAL  I  ZED  RESIDUAL 

N.B.  T  «  TREND  COMPONENT 

N.B.  E  *   OE- TRENDED,  DE-SESEASONALIZED  RESIDUAL 


SSS     VARIATION  IN  DE-TRENOED  RESIDUALS  SSS 

VARIATION  DUE  TO  SEASONAL  *  27.61154  X 
VARIATION  DUE  TO  NOISE  =  72.38847  X 
APPROXIMATE  RSO  ■  27.61153  X 


MEAN  X 
STD  DEV  X 


51.791650 
18.676070 


t«AN  E 
STD  DEV  E 


.310286 
16.134380 


ROBUST  SUMMARY  STATISTICS: 


MEDIAN  X    ■ 

51.109500 

MEDIAN  E 

.181818 

MEDIAN  |Xi   • 

51.109500 

MEDIAN  iE| 

7.465222 

MAO  X 

9.«?V94 

MAD  E 

7.283405 

.5  X  IM 

9.879261 

.5  E  IQR 

7.597422 

R  STD  DEV 

14.564330 

R  STO  DEV 

10.798390 

N.B.  X  >  DE-TRENOED  RESIDUAL 

N.B.  S  ■  SEASONAL  COMPONENT 

N.B.  E  '   DE-TRENDED,  DE-SESEASONALIZED  RESIDUAL 


APfQDIX:  lEPORT  FILES:  Tt2  mtmforMd  dita 
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SSS  TOTAL  VARIATION  IN  SERIES:  US 


VARIATION  DUE  TO  TREND 
VARIATION  DUE  TO  SEASON 
VARIATION  DUE  TO  S+T 
VARIATION  DUE  TO  NOISE 
APPROXIHATE  RSO 


59.83679  X 
-7.66329  X 
71.U390  X 
28.85609  X 
71.U391  X 


NEAN  X 
STD  DEV  X 


55.469480 
30.035590 


NEAN  E 
STD  DEV  E 


-.058188 
16.134380 


ROBUST   SmWARY   STATISTICS 


MEDIAN  X 
MEDIAN  jXJ 
MAD  X 
.5  X  lOR 
R  STD  DEV 


51.000000 
51.000000 
21.000000 
20.812500 
31.134650 


MEDIAN  E 
MEDIAN  JEJ 
MAC  E 
.5  E  IQR 
R  STD  DEV 


-.186657 

7.436587 

7.283407 

7.597423 

10.798400 


N.B.  X  >  RAW  DATA 

N.B.  T  X  TREND  COMPONENT 

N.B.  S  *   SEASONAL  COMPONENT 

N.B.  E  :  DE-TRENDED,  DE-SESEASONALIZED  RESIDUAL 

N.B.  D  =  CHRONOLOGICAL  SAMPLE  DATE 


SUH4ARY  STATISTICS  F0«  DE - SEASONAL  I  ZED ,  DE-TRENDED  RESIDUALS 
OeS  213 


A-MEAN 

-.0582 

STD  DEV 

16.1344 

MAX 

78.9488 

090 

16.8818 

075 

7.0805 

050 

-.1867 

Q25 

-8.1143 

010 

-17.0278 

MIN 

-80.1881 

RANGE  STATS 

MAX-MIN 

159.1370 

090-Q10 

33.9096 

075 -025 

15.1948 

TUKEY  FENCES  (STANDARDIZED)  FOR  OUTLIER  SCAN 


UPPER  OUTER  FENCE 

UPPER  INNER  FENCE 

LOWER  INNER  FENCE 

LOWER  OUTER  FENCE 


52.6651 

29.8728 

•30.9066 

-53.6988 


7.0819 
4.0170 
-4.1560 
-7.2209 


III  OUTLIERS  IN  DE-TRENOED  t  DE- SEASONAL  I  ZED  RESIDUALS  III 
DATE  TINE  RAW  DATUM   RESIDUAL  STANDARDIZED 


RESIDUAL  730809  1200 
RESIDUAL  741114  1200 
RESIDUAL  751014  1200 


162.000 

60.176 

8.092 

IS  FAR  OUT 

151.000 

78.949 

10.616 

IS  FAR  OUT 

91.000 

30.592 

4.114 

IS  OUT 

APf>EIBIX:  KPOIT  FILES:  TI2  XMttrwmformad  dita 
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RESIDUAL  760511  1200 
RESIDUAL  760607  1200 
RESIDUAL  780228  1200 
RESIDUAL  790710  930 
RESIDUAL  800730  1055 
RESIDUAL  870512  1200 
RESIDUAL  880607  1027 
RESIDUAL  890516  900 
RESIDUAL  890613  850 


126.000 

59.248 

7.967 

IS  FAR  OUT 

19.000 

-40.061 

-5.387 

IS  OUT 

20.000 

•80.188 

-10.783 

IS  FAR  OUT 

35.000 

-34.367 

-4.621 

IS  OUT 

U.OOO 

-49.704 

-6.684 

IS  ALMOST  FAR 

76.000 

37.768 

5.079 

IS  OUT 

50.000 

28.153 

3.786 

IS  ALMOST  OUT 

57.000 

29.195 

3.926 

IS  ALMOST  OUT 

51.000 

30.421 

4.091 

IS  OUT 

THERE  ARE:     2  ALMOST  OUT  DATA 
5  OUT  DATA 

1  ALMOST  FAR  OUT  DATA 
4  FAR  OUT  DATA 
12  TOTAL  OUTLIERS  IDENTIFIED 


N.B.   OUTLIERS  ARE  STANDARDIZED  BY  DIVISION  BY  THE  PCDIAN  ABSOLUTE  RESIDUAL  =  S 

FOR  HCRMAL  DISTRIBUTIONS: 
OUTER  FENCES  ARE  ABOUT  7S  OR  4.7  X  STANDARD  DEVIATION 
INNER  FENCES  ARE  ABOUT  4S  OR  2.7  X  STANDARD  DEVIATION 
ALSO  6S  OR  4  X  STANDARD  DEVIATION  (ALMOST  FAR  OUT) 
t   3.75S  OR  2.5  X  STANDARD  DEVIATION  (ALMOST  OUT) 
ARE  STANDARDIZED  VALUES  USE  BY  OTHER  AUTHORS 


TUICEY  FENCES  ARE  DEFINED  USING  1/2  THE  INTERQUARTILE  RAN(^ 
OF  RESIDUALS  WHICH  IS  NOT  EXACTLY  IDENTICAL  TO  S 
EXCEPT  FOR  NORMAL  DISTRIBUTIONS  BUT  IS  ALSO  A  COMMON 
ROeuST  MEASURE  OF  DISPERSION 


APPQOIX:  KPORT  FIL£S:  TI2  lO0^g  tnraforMd  dita  23 


PROGRAM  TR2 

REPORT  DATA  FILE  ^  D:\DlLDO\TR2\CRno3L 
INPUT  DATA  FILE  >  O:\DILD0\DAT2\N03.dat 
VARIABLE  -  M03 
(5X,i7,i5,f9.3) 
START  YEAR  >  720000   END  YEAR  -  920000 

ZERO  OPTION  CODE  -  1 

NUMERIC  MISSING  VALUE  CODE  <    -9.00000000 

MULTIPLIER/DIVISOR  CODE  -  1 

MULTIPLIER/DIVISOR  >      .00000000 

INPUT  DATA  ARE  LEFT  AS  IS 

LOGG  CODE  =   1 

RUNS  <   1 

ITERATIONS  >   5 

FILTER  «  5 

SAMPLING  INTERVAL  CODE  =  1 

MINIMUM  SMOOTHING  BASED  ON  MEDIAN  SAMPLING  INTERVAL 

UNIT  SMOOTHING  INTERVAL  =   1.00  yr 

MINIMUM  N  /  WINDOW  '      13 

MINIMUM  SEASONAL  SMOOTHING  WINDOW  -  1/  8  yr 

PLOT  COOES: 
PLOT  INDICATOR  >  0 
PLOT  DENSITY  1  =  D 
LOG  PLOT  TYPE  «  0 

PLOT  TYPE  CODE  «  L   (i.e.  LONG/SHORT) 
PLOT  HEADER  1  *   Nitrate  ug/L 
PLOT  HEADER  2  >  CLEARWATER  L 
RUN  DATE  HEADER  *   JUN  13  1991 
Y-AXIS  TITLE  '   LoglO  Concentration  ug/L 
REPORT  FILE  CODE  >  SR  (i.e.  SHORT  REPORT) 
OUTPUT  FILE  CODE  «  0 
HO  OUTPUT  FILES  REQUESTED 
FIRST  LINE  OF  INPUT  DATA  FILE: 

PARAMETER  «  N03      SITE  =  CLEARWATER  L 

QUALITY  DATA  FOR  PARM  ■  N03     READ  (  STORED 
213  RECORDS  WERE  READ 

THE  FOLLOWING  3  INPUT  DATA  WERE  IDENITIFIED  AS  ZERO'S: 
870901  1100     .000  ZERO 
880830  1030     .000  ZERO 
880926  1000      .000  ZERO 

SMALLEST  NON-ZERO  POSITIVE  VALUE  >     4.00000 

ZERO  OPTION  CODE  >  1 

SETTING  O'S  TO  1/2  SMALLEST  POSITIVE  NON-ZERO  VALUE  =     2.00000 


DOING  LOG10  TRANSFORM  OF  DATA  SERIES 


*««««««•« 


BEGINNING  RUN  1 


UNIT  SMOOTHING  INTERVAL  (USD  >      1.00000  yr 


LENGTH  OF  RECORD 
MEAN  SAMPLING  INTERVAL 
MEAN  SAMPLING  RATE/US  I 
MEDIAN  SAMPLING  INTERVAL 
MEDIAN  SAMPLING  RATE/USI 
MINIMUM  SMOOTHING 


16.60808  yr 
.07834  yr 

12.76487  #/USI 
.07382  yr 

13.54609  «/USI 
14  OBSERVATIONS 


APPQDIX:  lEPGKT  FILES:  TK2  I 


••10 


trwwforaad  data 
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TREND  SUMMARY: 
GRAND  TINE-UEICHTED  MEAN  TREND 


1.691 
49.068 


LOG  UNITS 
NATURAL  UNITS 


SCAN  SMOOTH  FOR  LOCAL  EXTREMA  INCLUDING  END  POINTS 


LOCAL  MAX  ' 
LOCAL  MIN  ■ 
LOCAL  MAX  ■ 
LOCAL  MIN  ' 
LOCAL  MAX  ■ 
LOCAL  MIN  ' 
LOCAL  MAX  c 
LOCAL  MIN  > 
LOCAL  MAX  * 
LOCAL  MIN  * 
LOCAL  MAX  * 
LOCAL  MIN  s 

MAX  TREND  > 
MIN  TREND  = 
TREND  RANGE 


148.124  8 
81 .997  a 
8S.250  a 
58.727  a 
94.778  a 
73.930  a 
77.887  a 
22.945  a 
29.177  a 
28.706  a 
37.119  a 
12.201  a 


73.489 
74.522 
74.870 
76.818 
78.656 
79.637 
80.788 
84.657 
85.483 
85.807 
86.960 
90.097 


JUN  28 
JUL  10 
NOV  14 
OCT  26 
AUG  28 
AUG  21 
OCT  15 
AUG  28 
JUN  26 
OCT  22 
DEC  17 
FEB  5 


1973 
1974 
1974 
1976 
1978 
1979 
1980 
1984 
1985 
1985 
1986 
1990 


148.124  a  73.489  >  JUN  28  1973 
12.201  a  90.097  >  FEB  5  1990 
135.923 


ANNUAL  TIME-WEIGHTED  MEAN  TREND  SUMMARY  ON 

LOG10  SER 

lES 

YEAR 

N 

ITER  «  1 

ITER  »  5 

DIFF 

YR-TO-YR 

DIFF 

1 

73 

16 

2.027 

2.131 

-.105 

74-73 

-.213 

2 

74 

21 

1.833 

1.918 

-.085 

75-74 

-.040 

3 

75 

10 

1.828 

1.879 

-.050 

76-75 

-.079 

4 

76 

13 

1.764 

1.800 

-.036 

77-76 

.035 

5 

77 

10 

1.809 

1.835 

-.026 

78-77 

.121 

6 

78 

14 

1.928 

1.956 

-.028 

79-78 

-.068 

7 

79 

19 

1.878 

1.888 

-.010 

80-79 

-.001 

8 

80 

12 

1.864 

1.887 

-.023 

81-80 

-.018 

9 

81 

7 

1.822 

1.869 

-.047 

82-81 

-.122 

10 

82 

7 

1.734 

1.746 

-.012 

83-82 

-.192 

11 

83 

10 

1.557 

1.555 

.002 

84-83 

-.169 

12 

84 

13 

1.385 

1.386 

-.001 

85-84 

.068 

13 

85 

12 

1.457 

1.454 

.003 

86-85 

.071 

14 

86 

12 

1.528 

1.525 

.004 

87-86 

-.016 

15 

87 

13 

1.543 

1.509 

.oil, 

88-87 

-.118 

16 

88 

13 

1.U7 

1.391 

.056 

89-88 

-.111 

17 

89 

10 

1.306 

1.280 

.026 

90-89 

-.194 

18 

90 

1 

1.059 

1.086 

-.027 

MEAN 

1  TREND 

- 

1.676 

1.672 

MAXIMUM  TREND  ' 

2.131  a 

1973 

MINIMUM  TREND  ■ 

1.086  a 

1990 

TREND  RANGE   ■ 

1.045 

ANNUAL  T-W  MEAN  TREND  SUMMARY  ON  NATURAL  SERIES 


YEAR 

N 

ITER  ■  1 

ITER  >  5 

DIFF 

YR-TO-YR 

DIFF 

73 

16 

106.395 

135.351 

-28.956 

74-73 

-52.527 

74 

21 

68.029 

82.824 

-14.795 

75-74 

-7.203 

75 

10 

67.335 

75.621 

-8.286 

76-75 

-12.511 

76 

13 

58.025 

63.110 

-5.085 

77-76 

5.276 

77 

10 

64.452 

68.385 

-3.934 

78-77 

21.998 

78 

14 

84.698 

90.384 

-5.686 

79-78 

-13.103 

79 

19 

75.537 

77.281 

-1.7U 

80-79 

-.202 

APKmit:   lEPatT  FILES:  TI2  log^g  trmforwd  diU 


8    80 

12 

73.169 

77.079 

-3.910 

81-80 

-3.150 

9    81 

7 

66.402 

73.929 

-7.527 

82-81 

-18.160 

10    82 

7 

54.260 

55.769 

-1.509 

83-82 

-19.905 

11    83 

10 

36.022 

35.864 

.158 

84-83 

-11.537 

12    84 

13 

24.247 

24.328 

-.081 

85-84 

4.115 

13    85 

12 

28.629 

28.443 

.187 

86-85 

5.028 

U    86 

12 

33.754 

33.471 

.283 

87-86 

-1.209 

15    87 

13 

34.926 

32.262 

2.664 

88-87 

-7.649 

16    88 

13 

28.022 

24.613 

3.409 

89-88 

-5.538 

17    89 

10 

20.239 

19.075 

1.164 

90-89 

-6.874 

18    90 

1 

11.456 

12.201 

-.7U 

G-MEAN  TREND  ' 

47.409 

49.068 

NAXIMUN 

TREND  ■ 

135.351  a 

1973 

MINIMUM 

TREND  « 

12.201  a 

1990 

TREND  RANGE   « 

123.150 

SPEARMANS  RHO  TREND  TEST  ON  FINAL  ANNUAL  TRENDS 

FOR  YEARS  WITH  WITH  5  OR  MORE  OBSERVATIONS 

IF  RHO  <  0  THEN  MET  TREND  IS  DOWN 

IF  RHO  >  0  THEN  MET  TREND  IS  UP 

RHO  '     -.86765   APPROXX  p  *      .000006   dof  *     IS 


SSSSSSS   SEASONAL  SUMMARY  SSSSSS 


MEAN  ANNUAL  SEASONAL  MEDIAN  ■ 
G-MEAN  ANNUAL  SEASONAL  MEDIAN  > 


1.690376  LOG  UNITS 
49.020320  NATURAL  UNITS 


SCAN  FOR  SEASONAL  EXTREMA 


LOCAL  MAX  < 
LOCAL  MIN  > 


.190614  a  .25  '     APR  2 
.131510  a  .62  >  AUG  15 


SEASONAL  MAX  >  .190614  a  .25  >  APR  2 
SEASONAL  MIN  >  -.131510  a  .62  «  AUG  15 
SEASONAL  AMPLITUDE  ■     .322124 


IN  NATURAL  UNITS 

AT  G-MEAN  ANNUAL  SEASONAL  MEDIAN  > 


49.020 


SEASONAL  MAX  «  76.030840  a  .25  ■  APR  2 
SEASONAL  MIN  «  36.213120  8  .62  «  AUG  15 
SEASONAL  AMPLITUDE  *       39.817720 


AT  MAXIMUM  TREND  «    135.351   IN  1973 

SEASONAL  MAX  >  209.930000  8  .25  ■  APR  2 
SEASONAL  MIN  >  99.988640  a  .62  ■  AUG  15 
SEASONAL  AMPLITUDE  ■  109.941400 


AT  MINIMUM  TREND  > 


12.201   IN  1990 


SEASONAL  MAX  «  18.923100  8  .25  «  APR  2 
SEASONAL  MIN  «  9.012978  8  .62  ■  AUG  15 
SEASONAL  AMPLITUDE  «    9.910118 


APPEIBIX:  IEI«T  FILES:  TU  Ioo^q  trwpmiormmi  datm 
ITERATIONS  >  5  COMPLETED 

SS$  VARIATION  IN  OE- SEASONAL  I  ZED  RESIDUALS  SSS 


VARIATION  DUE  TO  TREND 

62.19750  X 

VARIATION  DUE  TO  NOISE    ■ 

37.80250  X 

APPROXIMATE  RSQ          ■ 

62.19750  X 

MEAN  X      '            1.691450 

MEAN  E     ' 

-.029336 

STD  OEV  X   =      .326927 

STD  DEV  E   < 

.198M5 

ROBUST  SUMURT  STATISTICS: 


MEDIAN  X 

1.749538 

MEDIAN  E    > 

.oona 

MEDIAN  |X| 

1.749538 

MEDIAN  {E|  > 

.077308 

MAO  X 

.191625 

MAD  E 

.079055 

.5  X  lOR 

.205832 

.5  E  lOR 

.079137 

R  STD  DEV 

.284104 

R  STD  DEV   > 

.117207 

N.B.  X  -   DE-SESEASONALIZED  RESIDUAL 

N.B.  T  *   TREND  COMPONENT 

N.B.  E  <  DE-TRENOEO,  DE-SESEASONALIZED  RESIDUAL 


%SS     VARIATION  IN  OE-TRENOED  RESIDUALS  SSS 


VARIATION  DUE  TO  SEASONAL 

'     29.24598  X 

VARIATION  DUE  TO  NOISE 

■  70.75402  X 

APPROXIMATE  RSfi 

■  29.24598  X 

MEAN  X     «     1.625531 

MEAN  E      - 

-.028916 

STD  DEV  X   «     .229881 

STD  OEV  E  ' 

.198845 

ROBUST  SUMMARY  STATISTICS: 


MEDIAN  X 

1.641399 

MEDIAN  E 

X 

.003883 

MEDIAN  |xi 

1.641399 

MEDIAN  iEi 

S 

.077601 

MAD  X 

.104543 

MAD  E 

■ 

.079055 

.5  X  lOR 

.104333 

.5  E  lOR 

« 

.079137 

R  STD  DEV 

.154996 

R  STD  DEV 

s 

.117207 

N.B.  X  -  OE-TRENOED  RESIDUAL 

N.B.  S  ■  SEASONAL  COMPONENT 

N.B.  E  ■  DE-TRENOED,  DE-SESEASONALIZED  RESIDUAL 


SSS  TOTAL  VARIATION  IN  SERIES:  SSS 


VARIATION  DUE  TO  TREND 
VARIATION  DUE  TO  SEASON 

VARIATION  DUE  TO  S*T 
VARIATION  DUE  TO  NOISE 
APPROXIMATE  RSQ 


47.00081  X 
.78965  X 
62.50757  X 
37.49243  X 
62.50757  X 


APKIDIX:  KPOIT  FILES:  TI2  Ioq^q  tmforwd  dsta 
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MEAN  X 

■ 

1.655521 

MEAN  E 

■ 

-.029126 

STD  DEV  X 

« 

.326339 

STO  DEV  E 

■ 

.198045 

ROBUST  SUMMARY  STATISTICS 


POIAN  X 

1.707570 

MEDIAN  E 

.003673 

MEDIAN  |X; 

1.707570 

MEDIAN  JE|  • 

.077391 

MAD  X 

.167491 

MAD  E 

.079055 

.5  X  lOR 

.175905 

.5  E  lOR 

.079137 

R  STD  DEV 

.248323 

R  STD  DEV   > 

.117207 

N.B.  X  >  RAW  DATA 

N.B.  T  -  TREND  COMPONENT 

N.B.  S  ■  SEASONAL  COMPONENT 

N.B.  E  >  DE-TRENDED,  DE-SESEASONALIZEO  RESIDUAL 

N.B.  D  ■  CHRONOLOGICAL  SAMPLE  DATE 


SUMMARY  STATISTICS  FOR  DE- SEASONAL  I ZED,  DE- TRENDED  RESIDUALS 


213 


A-MEAN 

-.0291 

STD  DEV 

.1988 

MAX 

.4433 

090 

.1354 

075 

.0676 

050 

.0037 

025 

-.0907 

010 

-.1682 

MIN 

-1.0545 

RANGE  STATS 

MAX-MIN 

1.4978 

090-010 

.3035 

075-025 

.1583 

TUKEY  FENCES  (STANDARDIZED)  FOR  OUTLIER  SCAN 


UPPER 

OUTER  FENCE  ■ 

.5424 

7.0084 

UPPER 

INNER  FENCE  ■ 

.3050 

3.9407 

LOUER 

INNER.  FENCE  ■ 

-.3281 

-4.2398 

LOWER 

OUTER  FENCE  ■ 

-.5655 

-7.3075 

III  OUTLIERS  IN  DE-TREHOED  t  DE- SEASONAL  I  ZED  RESIDUALS  III 


DATE  TINE  RAW  DATUM 


LOG10   RESIDUAL  STANDARDIZED 


RESIDUAL  760607 

1200 

19.000 

1.279 

-.494 

-6.383 

IS 

ALMOST  FAR 

OUT 

RESIDUAL  780228 

1200 

20.000 

1.301 

-.760 

-9.815 

IS 

FAR  OUT 

RESIDUAL  800730 

1055 

14.000 

1.146 

-.624 

-8.065 

IS 

FAR  OUT 

RESIDUAL  821018 

1130 

20.000 

1.301 

-.316 

-4.083 

IS 

ALMOST  OUT 

RESIDUAL  840828 

955 

5.000 

.699 

-.535 

-6.908 

IS 

ALMOST  FAR 

OUT 

RESIDUAL  870512 

1200 

76.000 

1.881 

.325 

4.195 

IS 

OUT 

RESIDUAL  870901 

1100 

2.000 

.301 

-1.055 

-13.626 

IS 

FAR  OUT 

RESIDUAL  870928 

1010 

4.000 

.602 

-.765 

-9.885 

IS 

FAR  OUT 

RESIDUAL  880607 

1027 

50.000 

1.699 

.349 

4.508 

IS 

OUT 

RESIDUAL  880803 

1020 

5.000 

.699 

-.556 

-7.181 

IS 

ALMOST  FAR 

OUT 

RESIDUAL  880830 

1030 

2.000 

.301 

-.949 

-12.258 

IS 

FAR  OUT 

RESIDUAL  880926 

1000 

2.000 

.301 

-.968 

-12.509 

IS 

FAR  OUT 

RESIDUAL  890516 

900 

57.000 

1.756 

.414 

5.353 

IS 

OUT 

APPBOIX:  KPOtT  FILES:  TK2  Io0^q  tmforwd  dita  28 


RESIDUAL  890613  850     51.000 

1.708 

.U3 

5.727 

IS  OUT 

RESIDUAL  890711  830     41.000 

1.613 

.400 

5.167 

IS  OUT 

RESIDUAL  890925  1055     4.000 

.602 

-.509 

-6.573 

IS  ALMOST  FAR  CUT 

THERE  ARE:     1  ALMOST  OUT  DATA 

5  OUT  DATA 

4  ALMOST  FAR  CUT  DATA 

6  FAR  OUT  DATA 

16  TOTAL  OUTLIERS  IDENTIFIED 


N.B.  OUTLIERS  ARE  STANDARDIZED  BY  DIVISION  BY  THE  MEDIAN  ABSOLUTE  RESIDUAL  >  S 

FOR  NORMAL  DISTRIBUTIONS: 
OUTER  FENCES  ARE  ABOUT  7S  OR  4.7  X  STANDARD  DEVIATION 
INNER  FENCES  ARE  ABOUT  4S  OR  2.7  X  STANDARD  DEVIATION 
ALSO  6S  OR  4  X  STANDARD  DEVIATION  (ALMOST  FAR  OUT) 
t  3.75S  OR  2.5  X  STANDARD  DEVIATION  (ALMOST  OUT) 
ARE  STANDARDIZED  VALUES  USE  BY  OTHER  AUTHORS 

N.B.   TUKEY  FENCES  ARE  DEFINED  USING  1/2  THE  INTERQUARTILE  RANGE 
OF  RESIDUALS  WHICH  IS  NOT  EXACTLY  IDENTICAL  TO  S 
EXCEPT  FOR  NORMAL  DISTRIBUTIONS  BUT  IS  ALSO  A 
ROBUST  MEASURE  OF  DISPERSION 


